fc90b1992a06e2d510ff5a9d149cfd18d1d570f1
[u/mrichter/AliRoot.git] / HLT / src / AliL3Evaluate.cxx
1 //$Id$
2
3 // Author: Anders Vestbo <mailto:vestbo@fi.uib.no>
4 //*-- Copyright &copy ASV 
5
6 #include "AliL3StandardIncludes.h"
7 #include <TFile.h>
8 #include <TH1.h>
9 #include <TParticle.h>
10 #include <TTree.h>
11 #include <TClonesArray.h>
12 #include <TNtupleD.h>
13
14 #include <AliRun.h>
15 #include <AliSimDigits.h>
16 #include <AliTPC.h>
17 #include <AliTPCcluster.h>
18 #include <AliTPCClustersArray.h>
19 #include <AliTPCClustersRow.h>
20 #include <AliTPCParam.h>
21 #include <AliComplexCluster.h>
22
23 #if GCCVERSION == 3
24 #include <fstream>
25 #include <iosfwd>
26 #else
27 #include <fstream.h>
28 #endif
29
30 #include "AliL3Logging.h"
31 #include "AliL3Transform.h"
32 #include "AliL3SpacePointData.h"
33 #include "AliL3Track.h"
34 #include "AliL3FileHandler.h"
35 #include "AliL3TrackArray.h"
36 #include "AliL3Evaluate.h"
37
38 #if GCCVERSION == 3
39 using namespace std;
40 #endif
41
42 //_____________________________________________________________
43 // AliL3Evaluate
44 //
45 // Evaluation class for tracking. Plots efficiencies etc..
46 //
47
48
49 ClassImp(AliL3Evaluate)
50
51 AliL3Evaluate::AliL3Evaluate()
52 {
53   fTracks = 0;
54   fNFastPoints = 0;
55   fMcIndex = 0;
56   fMcId = 0;
57   fMinSlice=0;
58   fMaxSlice=0;
59   fGoodTracks=0;
60 }
61
62 AliL3Evaluate::AliL3Evaluate(Char_t *datapath,Int_t min_clusters,Int_t minhits,Double_t minpt,Int_t *slice)
63 {
64   
65   if(slice)
66     {
67       fMinSlice=slice[0];
68       fMaxSlice=slice[1];
69     }
70   else
71     {
72       fMinSlice=0;
73       fMaxSlice=35;
74     }
75   fMaxFalseClusters = 0.1;
76   fGoodFound = 0;
77   fGoodGen = 0;
78   fMinPointsOnTrack = min_clusters;
79   fMinHitsFromParticle = minhits;
80   fMinGoodPt = minpt;
81   Char_t fname[1024];
82   AliL3FileHandler *clusterfile[36][6];
83   for(Int_t s=fMinSlice; s<=fMaxSlice; s++)
84     {
85       for(Int_t p=0; p<AliL3Transform::GetNPatches(); p++)
86         {
87           fClusters[s][p] = 0;
88           clusterfile[s][p] = new AliL3FileHandler();
89           sprintf(fname,"%s/points_%d_%d.raw",datapath,s,p);
90           if(!clusterfile[s][p]->SetBinaryInput(fname))
91             {
92               LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
93                 <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; 
94               delete clusterfile[s][p];
95               clusterfile[s][p] = 0; 
96               continue;
97             }
98           fClusters[s][p] = (AliL3SpacePointData*)clusterfile[s][p]->Allocate();
99           clusterfile[s][p]->Binary2Memory(fNcl[s][p],fClusters[s][p]);
100           clusterfile[s][p]->CloseBinaryInput();
101         }
102     }
103
104   sprintf(fname,"%s/tracks.raw",datapath);
105   AliL3FileHandler *tfile = new AliL3FileHandler();
106   if(!tfile->SetBinaryInput(fname)){
107     LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
108       <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; 
109     return;
110   }
111   fTracks = new AliL3TrackArray();
112   tfile->Binary2TrackArray(fTracks);
113   tfile->CloseBinaryInput();
114   delete tfile;
115 }
116
117
118 AliL3Evaluate::~AliL3Evaluate()
119 {
120   if(fGoodTracks) delete fGoodTracks;
121   if(fDigitsTree) fDigitsTree->Delete();
122   if(fTracks) delete fTracks;
123   if(fPtRes) delete fPtRes;
124   if(fNGoodTracksPt) delete fNGoodTracksPt;
125   if(fNFoundTracksPt) delete fNFoundTracksPt;
126   if(fNFakeTracksPt) delete fNFakeTracksPt;
127   if(fTrackEffPt) delete fTrackEffPt;
128   if(fFakeTrackEffPt) delete fFakeTrackEffPt;
129   if(fNGoodTracksEta) delete fNGoodTracksEta;
130   if(fNFoundTracksEta) delete fNFoundTracksEta;
131   if(fNFakeTracksEta) delete fNFakeTracksEta;
132   if(fTrackEffEta) delete fTrackEffEta;
133   if(fFakeTrackEffEta) delete fFakeTrackEffEta;
134   if(fMcIndex) delete [] fMcIndex;
135   if(fMcId)    delete [] fMcId;
136   if(fNtuppel) delete fNtuppel;
137 }
138
139
140 void AliL3Evaluate::AssignIDs()
141 {
142   //Assign MC id to the tracks.
143 #ifndef do_mc
144   cerr<<"AliL3Evaluate::AssignIDs() : You need to compile with the do_mc flag!"<<endl;
145   return;
146 #endif
147   
148   fTracks->QSort();
149   LOG(AliL3Log::kDebug,"AliL3Evaluate::AssignIDs","Track Loop")
150     <<"Assigning MC id to the found tracks...."<<ENDLOG;
151   for(Int_t i=0; i<fTracks->GetNTracks(); i++)
152     {
153       AliL3Track *track = (AliL3Track*)fTracks->GetCheckedTrack(i);
154       if(!track) continue; 
155       if(track->GetNumberOfPoints() < fMinPointsOnTrack) break;
156       
157       fGoodFound++;
158       Int_t tID = GetMCTrackLabel(track);
159       track->SetMCid(tID);
160     }
161   cout<<"Found "<<fGoodFound<<" good tracks "<<endl;
162 }
163
164
165 struct S {Int_t lab; Int_t max;};
166 Int_t AliL3Evaluate::GetMCTrackLabel(AliL3Track *track){ 
167   //Returns the MCtrackID of the belonging clusters.
168   //If MCLabel < 0, means that track is fake.
169   //Fake track means that more than 10 percent of clusters are assigned incorrectly.
170   
171 #ifdef do_mc
172   Int_t num_of_clusters = track->GetNumberOfPoints();
173   S *s=new S[num_of_clusters];
174   Int_t i;
175   for (i=0; i<num_of_clusters; i++) s[i].lab=s[i].max=0;
176   UInt_t *hitnum = track->GetHitNumbers();  
177   UInt_t id;
178     
179   Int_t lab=123456789;
180   for (i=0; i<num_of_clusters; i++) 
181     {
182       //Tricks to get the clusters belonging to this track:
183       id = hitnum[i];
184       Int_t slice = (id>>25) & 0x7f;
185       Int_t patch = (id>>22) & 0x7;
186       UInt_t pos = id&0x3fffff;       
187       
188       AliL3SpacePointData *points = fClusters[slice][patch];
189       if(!points) continue; 
190       if(pos>=fNcl[slice][patch]) 
191         {
192           LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray")
193             <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
194           continue;
195         }
196       
197       //Get the label of the cluster:
198       
199       lab=abs(points[pos].fTrackID[0]);
200       if(lab < 0)
201         cout<<"Track had negative id : "<<lab<<" padrow "<<(Int_t)points[pos].fPadRow<<" nhits "<<num_of_clusters<<" pt "<<track->GetPt()<<endl;
202       
203       Int_t j;
204       for (j=0; j<num_of_clusters; j++)
205         if (s[j].lab==lab || s[j].max==0) break;
206       s[j].lab=lab;
207       s[j].max++;
208     }
209   
210   Int_t max=0;
211   for (i=0; i<num_of_clusters; i++) 
212     if (s[i].max>max) {max=s[i].max; lab=s[i].lab;}
213   
214   delete[] s;
215   
216   for (i=0; i<num_of_clusters; i++) 
217     {
218       id = hitnum[i];
219       Int_t slice = (id>>25) & 0x7f;
220       Int_t patch = (id>>22) & 0x7;
221       UInt_t pos = id&0x3fffff;       
222       
223       AliL3SpacePointData *points = fClusters[slice][patch];
224       if(!points) continue; 
225       if(pos>=fNcl[slice][patch]) 
226         {
227           LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray")
228             <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
229           continue;
230         }
231       
232       if (abs(points[pos].fTrackID[1]) == lab || 
233           abs(points[pos].fTrackID[2]) == lab ) max++;
234     }
235   
236   if (1.-Float_t(max)/num_of_clusters > fMaxFalseClusters) 
237     {
238       return -lab;
239     }
240   
241   return lab;
242 #else //If we are running with mc_ids or not
243   return 0;
244 #endif
245
246 }
247
248 void AliL3Evaluate::GetFastClusterIDs(Char_t *path)
249 {
250   //Get the MC id of space points in case of using the fast simulator. 
251   char fname[256];
252   sprintf(fname,"%s/point_mc.dat",path);
253   FILE *infile = fopen(fname,"r");
254   if(!infile) return;
255   Int_t hitid,hitmc,i;
256   
257   for(i=0; ; i++)
258     if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break;
259   rewind(infile);
260   fNFastPoints = i;
261   fMcId = new Int_t[fNFastPoints];
262   fMcIndex = new UInt_t[fNFastPoints];
263   
264   for(i=0; i<fNFastPoints; i++)
265     {
266       if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break;
267       fMcId[i] = hitmc;
268       fMcIndex[i] = hitid;
269     }
270   fclose(infile);
271 }
272
273 void AliL3Evaluate::CreateHistos(Int_t nbin,Float_t xlow,Float_t xup)
274 {
275   //Create the histograms 
276   
277   LOG(AliL3Log::kInformational,"AliL3Evaluate::CreateHistos","Allocating")
278     <<"Creating histograms..."<<ENDLOG;
279   
280   fNtuppel = new TNtuple("fNtuppel","Pt resolution","pt_gen:pt_found:nHits");
281   fNtuppel->SetDirectory(0);
282   fPtRes = new TH1F("fPtRes","Relative Pt resolution",30,-10.,10.); 
283   fNGoodTracksPt = new TH1F("fNGoodTracksPt","Good tracks vs pt",nbin,xlow,xup);    
284   fNFoundTracksPt = new TH1F("fNFoundTracksPt","Found tracks vs pt",nbin,xlow,xup);
285   fNFakeTracksPt = new TH1F("fNFakeTracksPt","Fake tracks vs pt",nbin,xlow,xup);
286   fTrackEffPt = new TH1F("fTrackEffPt","Tracking efficiency vs pt",nbin,xlow,xup);
287   fFakeTrackEffPt = new TH1F("fFakeTrackEffPt","Efficiency for fake tracks vs pt",nbin,xlow,xup);
288   
289   fNGoodTracksEta = new TH1F("fNGoodTracksEta","Good tracks vs eta",20,-50,50);
290   fNFoundTracksEta = new TH1F("fNFoundTracksEta","Found tracks vs eta",20,-50,50);
291   fNFakeTracksEta = new TH1F("fNFakeTracksEta","Fake tracks vs eta",20,-50,50);
292   fTrackEffEta = new TH1F("fTrackEffEta","Tracking efficienct vs eta",20,-50,50);
293   fFakeTrackEffEta = new TH1F("fFakeTrackEffEta","Efficiency for fake tracks vs eta",20,-50,50);
294
295 }
296
297 void AliL3Evaluate::GetGoodParticles(Char_t *path,Bool_t sector)
298 {
299   //Read the good particles from file. This file should already have been
300   //generated by macro AliTPCComparison.C.
301   
302   Char_t filename[1024];
303   if(!sector)
304     sprintf(filename,"%s/good_tracks_tpc",path);
305   else
306     sprintf(filename,"%s/good_tracks_tpc_sector",path);//Sectorwise comparison.
307   ifstream in(filename);
308   if(!in)
309     {
310       cerr<<"AliL3Evaluate::GetGoodParticles : Problems opening file :"<<filename<<endl;
311       return;
312     }
313   Int_t MaxTracks=20000;
314   fGoodTracks = new GoodTrack[MaxTracks];
315   
316   while (in>>fGoodTracks[fGoodGen].label>>fGoodTracks[fGoodGen].code>>
317          fGoodTracks[fGoodGen].px>>fGoodTracks[fGoodGen].py>>fGoodTracks[fGoodGen].pz>>
318          fGoodTracks[fGoodGen].x>>fGoodTracks[fGoodGen].y >>fGoodTracks[fGoodGen].z>>fGoodTracks[fGoodGen].nhits>>fGoodTracks[fGoodGen].sector) 
319
320     {
321       fGoodGen++;
322       if (fGoodGen==MaxTracks) 
323         {
324           cerr<<"AliL3Evaluate::GetGoodParticles : Too many good tracks !\n";
325           break;
326         }
327     }  
328 }
329
330
331 void AliL3Evaluate::FillEffHistos()
332 {  
333   if(!fGoodTracks)
334     {
335       cerr<<"AliL3Evaluate::FillEffHistos : No good tracks"<<endl;
336       return;
337     }
338   cout<<"Comparing "<<fGoodGen<<" good tracks ..."<<endl;
339   for(Int_t i=0; i<fGoodGen; i++)
340     {
341       //cout<<"Checking particle "<<i<<endl;
342       if(fGoodTracks[i].nhits < fMinHitsFromParticle) continue;
343       Float_t ptg = TMath::Sqrt(fGoodTracks[i].px*fGoodTracks[i].px + fGoodTracks[i].py*fGoodTracks[i].py);
344       if(ptg < fMinGoodPt) continue;
345       Float_t pzg=fGoodTracks[i].pz;
346       Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi();
347
348       fNGoodTracksPt->Fill(ptg);
349       fNGoodTracksEta->Fill(dipangle);
350       Int_t found = 0;
351       
352       for(Int_t k=0; k<fTracks->GetNTracks(); k++)
353         {
354           AliL3Track *track = fTracks->GetCheckedTrack(k);
355           if(!track) continue;
356           Int_t nHits = track->GetNumberOfPoints();
357           if(nHits < fMinPointsOnTrack) break;
358           Int_t tracklabel;
359           tracklabel = track->GetMCid();
360           
361           if(TMath::Abs(tracklabel) != fGoodTracks[i].label) continue;
362           found=1;
363           if(tracklabel == fGoodTracks[i].label) {fNFoundTracksPt->Fill(ptg); fNFoundTracksEta->Fill(dipangle);}
364           else {fNFakeTracksPt->Fill(ptg); fNFakeTracksEta->Fill(dipangle);}
365           Float_t pt=track->GetPt();
366           fPtRes->Fill((pt-ptg)/ptg*100.);
367           fNtuppel->Fill(ptg,pt,nHits);
368           break;
369           
370         }
371       if(!found)
372         cout<<"Track "<<fGoodTracks[i].label<<" was not found"<<endl;
373     }
374 }
375
376 void AliL3Evaluate::FillEffHistosNAIVE()
377 {  
378   //Fill the efficiency histograms.
379   
380   cout<<endl<<"Note: Doing NAIVE evaluation "<<endl;
381   for(Int_t i=0; i<fGoodGen; i++)
382     {
383       Double_t ptg=TMath::Sqrt(fGoodTracks[i].px*fGoodTracks[i].px + fGoodTracks[i].py*fGoodTracks[i].py);
384       Double_t pzg=fGoodTracks[i].pz;
385       Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi();
386       //printf("filling particle with pt %f and dipangle %f\n",ptg,dipangle);
387       fNGoodTracksPt->Fill(ptg);
388       fNGoodTracksEta->Fill(dipangle);
389       
390     }
391   
392   for(Int_t k=0; k<fTracks->GetNTracks(); k++)
393     {
394       AliL3Track *track = fTracks->GetCheckedTrack(k);
395       if(!track) continue;
396       Int_t nHits = track->GetNumberOfPoints();
397       if(nHits < fMinPointsOnTrack) break;
398       if(track->GetPt()<fMinGoodPt) continue;
399       if(fabs(track->GetPseudoRapidity())>0.9) continue;
400
401       fNFoundTracksPt->Fill(track->GetPt()); fNFoundTracksEta->Fill(track->GetPseudoRapidity());
402       //Float_t pt=track->GetPt();
403       //fPtRes->Fill((pt-ptg)/ptg*100.);
404       //fNtuppel->Fill(ptg,pt,nHits);
405             
406     }
407 }
408
409 void AliL3Evaluate::CalcEffHistos(){  
410   
411   Stat_t ngood=fNGoodTracksPt->GetEntries();
412   Stat_t nfound=fNFoundTracksPt->GetEntries();
413   Stat_t nfake=fNFakeTracksPt->GetEntries();
414   
415   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
416     <<AliL3Log::kDec<<"There was "<<ngood<<" generated good tracks"<<ENDLOG;
417   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
418     <<AliL3Log::kDec<<"Found "<<nfound<<" tracks"<<ENDLOG;
419   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
420     <<AliL3Log::kDec<<"Integral efficiency is about "<<nfound/ngood*100<<ENDLOG;
421   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
422     <<AliL3Log::kDec<<"Fake tracks relative is about "<<nfake/ngood*100<<ENDLOG;
423   
424   LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency")
425     <<AliL3Log::kDec<<"Naive efficiency "<<(Double_t)fGoodFound/(Double_t)fGoodGen<<ENDLOG;
426
427   
428   fNFoundTracksPt->Sumw2(); fNGoodTracksPt->Sumw2();
429   fTrackEffPt->Divide(fNFoundTracksPt,fNGoodTracksPt,1,1.,"b");
430   fFakeTrackEffPt->Divide(fNFakeTracksPt,fNGoodTracksPt,1,1.,"b");
431   fTrackEffPt->SetMaximum(1.4);
432   fTrackEffPt->SetXTitle("P_{T} [GeV]");
433   fTrackEffPt->SetLineWidth(2);
434   fFakeTrackEffPt->SetFillStyle(3013);
435   fTrackEffPt->SetLineColor(4);
436   fFakeTrackEffPt->SetFillColor(2);
437   
438   fNFoundTracksEta->Sumw2(); fNGoodTracksEta->Sumw2();
439   fTrackEffEta->Divide(fNFoundTracksEta,fNGoodTracksEta,1,1.,"b");
440   fFakeTrackEffEta->Divide(fNFakeTracksEta,fNGoodTracksEta,1,1.,"b");
441   fTrackEffEta->SetMaximum(1.4);
442   fTrackEffEta->SetXTitle("#lambda [degrees]");
443   fTrackEffEta->SetLineWidth(2);
444   fFakeTrackEffEta->SetFillStyle(3013);
445   fTrackEffEta->SetLineColor(4);
446   fFakeTrackEffEta->SetFillColor(2);
447        
448 }
449
450 void AliL3Evaluate::Write2File(Char_t *outputfile)
451 {
452   //Write histograms to file:
453   
454   TFile *of = TFile::Open(outputfile,"RECREATE");
455   if(!of->IsOpen())
456     {
457       LOG(AliL3Log::kError,"AliL3Evaluate::Write2File","File Open")
458         <<"Problems opening rootfile"<<ENDLOG;
459       return;
460     }
461   
462   of->cd();
463   fNtuppel->Write();
464   fPtRes->Write();
465   fNGoodTracksPt->Write();
466   fNFoundTracksPt->Write();
467   fNFakeTracksPt->Write();
468   fTrackEffPt->Write();
469   fFakeTrackEffPt->Write();
470   fNGoodTracksEta->Write();
471   fNFoundTracksEta->Write();
472   fNFakeTracksEta->Write();
473   fTrackEffEta->Write();
474   fFakeTrackEffEta->Write();
475   
476   of->Close();
477
478 }
479
480 TNtupleD *AliL3Evaluate::CalculateResiduals(Char_t *datapath)
481 {
482
483   TNtupleD *ntuppel=new TNtupleD("ntuppel","Residuals","residual_trans:residual_long:zHit:pt:dipangle:beta:padrow:nHits");
484   ntuppel->SetDirectory(0);
485   
486   for(int f=fMinSlice; f<=fMaxSlice; f++)
487     {
488       AliL3FileHandler *tfile = new AliL3FileHandler();
489       char fname[256];
490       sprintf(fname,"%s/tracks_tr_%d_0.raw",datapath,f);
491       if(!tfile->SetBinaryInput(fname)){
492         LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open")
493           <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; 
494         return 0;
495       }
496       fTracks = new AliL3TrackArray();
497       tfile->Binary2TrackArray(fTracks);
498       tfile->CloseBinaryInput();
499       delete tfile;
500       printf("Looking in slice %d\n",f);
501       for(Int_t i=0; i<fTracks->GetNTracks(); i++)
502         {
503           
504           AliL3Track *track = (AliL3Track*)fTracks->GetCheckedTrack(i);
505           if(!track) continue;
506           if(track->GetNHits() < fMinPointsOnTrack) continue;
507           
508           track->CalculateHelix();
509           UInt_t *hitnum = track->GetHitNumbers();
510           UInt_t id;
511           
512           Float_t xyz[3];
513           Int_t padrow;
514           for(Int_t j=0; j<track->GetNumberOfPoints()-1; j++)
515             {
516               id = hitnum[j];
517               Int_t slice = (id>>25) & 0x7f;
518               Int_t patch = (id>>22) & 0x7;
519               UInt_t pos = id&0x3fffff;       
520               
521               //if(slice!=1) continue;
522               
523               AliL3SpacePointData *points = fClusters[slice][patch];
524               
525               if(!points) 
526                 {
527                   LOG(AliL3Log::kError,"AliL3Evaluate::CalculateResiduals","Clusterarray")
528                     <<"No points at slice "<<slice<<" patch "<<patch<<" pos "<<pos<<ENDLOG;
529                   continue;
530                 }
531               if(pos>=fNcl[slice][patch]) 
532                 {
533                   LOG(AliL3Log::kError,"AliL3Evaluate::CalculateResiduals","Clusterarray")
534                     <<AliL3Log::kDec<<"ERROR"<<ENDLOG;
535                   continue;
536                 }
537               
538               xyz[0] = points[pos].fX;
539               xyz[1] = points[pos].fY;
540               xyz[2] = points[pos].fZ;
541               padrow = points[pos].fPadRow;
542               AliL3Transform::Global2Local(xyz,slice);
543               
544               Float_t xyz_cross[3];
545               track->GetCrossingPoint(padrow,xyz_cross);
546               Double_t beta = track->GetCrossingAngle(padrow);
547               
548               Double_t yres = xyz_cross[1] - xyz[1];
549               Double_t zres = xyz_cross[2] - xyz[2];
550               
551               Double_t dipangle = atan(track->GetTgl());
552               ntuppel->Fill(yres,zres,xyz_cross[2],track->GetPt(),dipangle,beta,padrow,track->GetNumberOfPoints());
553               
554             }
555         }
556       if(fTracks)
557         delete fTracks;
558     }
559   return ntuppel;
560 }
561
562 TNtuple *AliL3Evaluate::EvaluatePoints(Char_t *rootfile)
563 {
564   //Compare points to the exact crossing points of track and padrows.
565   //The input file to this function, contains the exact clusters calculated
566   //in AliTPC::Hits2ExactClusters.
567     
568   cout<<"Evaluating points"<<endl;
569   TNtuple *ntuppel = new TNtuple("ntuppel","residuals","slice:padrow:resy:resz:zHit:pt");
570   
571   TFile *exfile = TFile::Open(rootfile);
572   if(!exfile)
573     {
574       cerr<<"Error opening rootfile "<<rootfile<<endl;
575       return 0;
576     }
577
578   AliTPCParam *param = (AliTPCParam*)exfile->Get(AliL3Transform::GetParamName());
579   
580   //Get the exact clusters from file:
581   AliTPCClustersArray *arr = new AliTPCClustersArray;
582   arr->Setup(param);
583   arr->SetClusterType("AliComplexCluster");
584   char treeName[500];
585   sprintf(treeName,"TreeCExact_%s",param->GetTitle());
586   Bool_t clusterok = arr->ConnectTree(treeName);//Segment Tree (for offline clusters)
587   if(!clusterok) {printf("AliL3Evaluate::EvaluatePoints : Error in clusterloading\n"); return 0;}
588   
589   cout<<"Entering loop with "<<(Int_t)arr->GetTree()->GetEntries()<<endl;
590   for(Int_t i=0; i<arr->GetTree()->GetEntries(); i++)
591     {
592       //Get the exact clusters for this row:
593       Int_t cursec,currow;
594       AliSegmentID *s = arr->LoadEntry(i);
595       param->AdjustSectorRow(s->GetID(),cursec,currow);
596       
597       AliTPCClustersRow *ro = (AliTPCClustersRow *)arr->GetRow(cursec,currow);
598       TClonesArray *clusters = ro->GetArray();
599       int num_of_offline=clusters->GetEntriesFast();
600       
601       //Get the found clusters:
602       Int_t slice,padrow;
603       AliL3Transform::Sector2Slice(slice,padrow,cursec,currow);
604       if(slice<fMinSlice || slice>fMaxSlice) continue;
605       AliL3SpacePointData *points = fClusters[slice][0];
606       
607       Int_t index = fRowid[slice][padrow];
608       if(!fDigitsTree->GetEvent(index))
609         printf("AliL3Evaluate::EvaluatePoints : ERROR IN DIGITSTREE\n");
610       printf("Checking slice %d padrow %d with %d clusters\n",slice,padrow,num_of_offline);
611       
612       for(UInt_t c=0; c<fNcl[slice][0]; c++)
613         {
614           if(points[c].fPadRow!=padrow) continue;
615           for(Int_t m=0; m<num_of_offline; m++)
616             {
617               AliComplexCluster *cluster = (AliComplexCluster *)clusters->UncheckedAt(m);
618               Int_t mcId = cluster->fTracks[0];
619               //Int_t mcId = cluster->GetLabel(0);
620               if(mcId <0) continue;
621               TParticle *part = gAlice->Particle(mcId);
622               
623               Float_t xyz_cl[3] = {points[c].fX,points[c].fY,points[c].fZ};
624               Float_t xyz_ex[3];
625               AliL3Transform::Global2Raw(xyz_cl,cursec,currow);
626               if(fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),0)!=mcId &&
627                  fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),1)!=mcId &&
628                  fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),2)!=mcId)
629                 continue;
630               AliL3Transform::Raw2Local(xyz_ex,cursec,currow,cluster->fY,cluster->fX);
631               AliL3Transform::Raw2Local(xyz_cl,cursec,currow,xyz_cl[1],xyz_cl[2]);
632               Float_t resy = xyz_cl[1] - xyz_ex[1];//cluster->GetY()
633               Float_t resz = xyz_cl[2] - xyz_ex[2];//cluster->GetZ()
634               
635               ntuppel->Fill(slice,padrow,resy,resz,xyz_ex[2],part->Pt());
636             }
637         }      
638       arr->ClearRow(cursec,currow);
639     }
640
641   return ntuppel;
642 }
643
644 void AliL3Evaluate::GetCFeff(Char_t *outfile)
645 {
646   
647   TNtuple *ntuppel = new TNtuple("ntuppel","Cluster finder efficiency","row:ncrossings:nclusters");
648   
649   AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC");
650   
651   TPC->SetParam(fParam);
652   
653   Int_t ver = TPC->IsVersion();
654   LOG(AliL3Log::kInformational,"AliL3Evaluate::GetCFeff","TPC version")
655     <<"TPC version "<<ver<<" found on file"<<ENDLOG;
656   
657   Int_t zero=TPC->GetParam()->GetZeroSup();
658   
659   Int_t np = gAlice->GetNtrack();
660   
661   Int_t crossed,recs;
662   Int_t *count = new Int_t[np]; //np number of particles.
663   Int_t i;
664   Float_t xyz[3];
665   for (i=0; i<np; i++) count[i]=0;
666   for(Int_t sl=fMinSlice; sl<=fMaxSlice; sl++)
667     {
668       for (i=0; i<=175; i++) 
669         {
670           crossed=0;
671           recs=0;
672           Int_t index = fRowid[sl][i];
673           if (!fDigitsTree->GetEvent(index)) continue;
674           Int_t sec,row;
675           fParam->AdjustSectorRow(fDigits->GetID(),sec,row);
676           fDigits->First();
677           do {
678             Int_t it=fDigits->CurrentRow(), ip=fDigits->CurrentColumn();
679             Short_t dig = fDigits->GetDigit(it,ip);
680             
681             if(dig<=fParam->GetZeroSup()) continue;
682             if(it < fParam->GetMaxTBin()-1 && it > 0)
683               if(fDigits->GetDigit(it+1,ip) <= fParam->GetZeroSup()
684                  && fDigits->GetDigit(it-1,ip) <= fParam->GetZeroSup())
685                 continue;
686             
687             AliL3Transform::Raw2Local(xyz,sec,row,ip,it);
688             if(fParam->GetPadRowRadii(sec,row)<230./250.*fabs(xyz[2]))
689               continue;
690             
691             
692             Int_t idx0=fDigits->GetTrackID(it,ip,0); 
693             Int_t idx1=fDigits->GetTrackID(it,ip,1);
694             Int_t idx2=fDigits->GetTrackID(it,ip,2);
695             
696             if (idx0>=0 && dig>=zero) count[idx0]+=1;
697             if (idx1>=0 && dig>=zero) count[idx1]+=1;
698             if (idx2>=0 && dig>=zero) count[idx2]+=1;
699           } while (fDigits->Next());
700           
701           for (Int_t j=0; j<np; j++) 
702             {
703               if (count[j]>1) //at least two digits at this padrow 
704                 {
705                   crossed++;
706                   count[j]=0;
707                 }
708             }
709           AliL3SpacePointData *points = fClusters[sl][0];
710           for(UInt_t k=0; k<fNcl[sl][0]; k++)
711             {
712               if(points[k].fPadRow!=i) continue;
713               recs++;
714             }
715           ntuppel->Fill(i,crossed,recs);
716         }
717       
718     }
719   delete[] count;
720   
721   TFile *file = TFile::Open(outfile,"RECREATE");
722   ntuppel->Write();
723   file->Close();
724   
725 }
726