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