Removing memory leaks
[u/mrichter/AliRoot.git] / TOF / AliTOFtracker.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15 // AliTOFtracker Class
16 // Task: Perform association of the ESD tracks to TOF Clusters
17 // and Update ESD track with associated TOF Cluster parameters 
18 //
19 // -- Authors : S. Arcelli, C. Zampolli (Bologna University and INFN) 
20 // -- Contacts: Annalisa.De.Caro@cern.ch
21 // --         : Chiara.Zampolli@bo.infn.it
22 // --         : Silvia.Arcelli@bo.infn.it
23 //--------------------------------------------------------------------
24
25 #include "AliTOFtracker.h"
26 #include "AliTOFtrack.h"
27 #include "TClonesArray.h"
28 #include "TError.h"
29 #include "AliTOFdigit.h"
30 #include "AliTOFGeometry.h"
31 #include "AliTOF.h"
32 #include "AliRun.h"
33 #include "AliModule.h"
34
35 ClassImp(AliTOFtracker)
36
37 //_____________________________________________________________________________
38 AliTOFtracker::AliTOFtracker(AliTOFGeometry * geom, Double_t parPID[2]) { 
39   //AliTOFtracker main Ctor
40
41   fHoles=true;
42   fNseeds=0;
43   fNseedsTOF=0;
44   fngoodmatch=0;
45   fnbadmatch=0;
46   fnunmatch=0;
47   fnmatch=0;
48   fGeom = geom;
49   fTOFpid = new AliTOFpidESD(parPID);
50   fR=378.; 
51   fTOFHeigth=15.3;  
52   fdCut=3.; 
53   fDy=AliTOFGeometry::XPad(); 
54   fDz=AliTOFGeometry::ZPad(); 
55   fDx=1.5; 
56   fSeeds=0x0;
57   fTracks=0x0;
58   fN=0;
59   Init(); // temporary solution to know about Holes/no Holes
60 }
61 //_____________________________________________________________________________
62 AliTOFtracker::AliTOFtracker(const AliTOFtracker &t):AliTracker() { 
63   //AliTOFtracker copy Ctor
64
65   fHoles=t.fHoles;
66   fNseeds=t.fNseeds;
67   fNseedsTOF=t.fNseedsTOF;
68   fngoodmatch=t.fngoodmatch;
69   fnbadmatch=t.fnbadmatch;
70   fnunmatch=t.fnunmatch;
71   fnmatch=t.fnmatch;
72   fGeom = t.fGeom;
73   fTOFpid = t.fTOFpid;
74   fR=t.fR; 
75   fTOFHeigth=t.fTOFHeigth;  
76   fdCut=t.fdCut; 
77   fDy=t.fDy; 
78   fDz=t.fDz; 
79   fDx=1.5; 
80   fSeeds=t.fSeeds;
81   fTracks=t.fTracks;
82   fN=t.fN;
83 }
84 //_____________________________________________________________________________
85 void AliTOFtracker::Init() { 
86
87 // temporary solution to know about Holes/no Holes, will be implemented as 
88 // an AliTOFGeometry getter
89
90   AliModule* frame=gAlice->GetModule("FRAME"); 
91
92   if(!frame) {
93     Error("Init","Could Not load FRAME! Assume Frame with Holes \n");
94     fHoles=true;
95   } else{
96     if(frame->IsVersion()==1) {fHoles=false;}    
97     else {fHoles=true;}      
98   }
99 }
100 //_____________________________________________________________________________
101 Int_t AliTOFtracker::PropagateBack(AliESD* event) {
102   //
103   // Gets seeds from ESD event and Match with TOF Clusters
104   //
105
106
107   //Initialise some counters
108
109   fNseeds=0;
110   fNseedsTOF=0;
111   fngoodmatch=0;
112   fnbadmatch=0;
113   fnunmatch=0;
114   fnmatch=0;
115
116   Int_t ntrk=event->GetNumberOfTracks();
117   fNseeds = ntrk;
118   fSeeds= new TClonesArray("AliESDtrack");
119   TClonesArray &aESDTrack = *fSeeds;
120
121
122   //Load ESD tracks into a local Array of ESD Seeds
123
124   for (Int_t i=0; i<fNseeds; i++) {
125     AliESDtrack *t=event->GetTrack(i);
126     new(aESDTrack[i]) AliESDtrack(*t);
127   }
128
129   //Prepare ESD tracks candidates for TOF Matching
130   CollectESD();
131
132   //First Step with Strict Matching Criterion
133   MatchTracks(kFALSE);
134
135   //Second Step with Looser Matching Criterion
136   MatchTracks(kTRUE);
137
138   Info("PropagateBack","Number of matched tracks: %d",fnmatch);
139   Info("PropagateBack","Number of good matched tracks: %d",fngoodmatch);
140   Info("PropagateBack","Number of bad  matched tracks: %d",fnbadmatch);
141
142   //Update the matched ESD tracks
143
144   for (Int_t i=0; i<ntrk; i++) {
145     AliESDtrack *t=event->GetTrack(i);
146     AliESDtrack *seed =(AliESDtrack*)fSeeds->UncheckedAt(i);
147     if(seed->GetTOFsignal()>0){
148       t->SetTOFsignal(seed->GetTOFsignal());
149       t->SetTOFcluster(seed->GetTOFcluster());
150       AliTOFtrack *track = new AliTOFtrack(*seed); 
151       t->UpdateTrackParams(track,AliESDtrack::kTOFout);    
152       delete track;
153     }
154   }
155
156
157   //Make TOF PID
158   fTOFpid->MakePID(event);
159
160   if (fSeeds) {
161     fSeeds->Delete();
162     delete fSeeds;
163     fSeeds = 0x0;
164   }
165   if (fTracks) {
166     fTracks->Delete();
167     delete fTracks;
168     fTracks = 0x0;
169   }
170   return 0;
171   
172 }
173 //_________________________________________________________________________
174 void AliTOFtracker::CollectESD() {
175    //prepare the set of ESD tracks to be matched to clusters in TOF
176  
177   fTracks= new TClonesArray("AliTOFtrack");
178   TClonesArray &aTOFTrack = *fTracks;
179   for (Int_t i=0; i<fNseeds; i++) {
180
181     AliESDtrack *t =(AliESDtrack*)fSeeds->UncheckedAt(i);
182     if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
183
184     // TRD good tracks, already propagated at 371 cm
185
186     if (
187         ((t->GetStatus()&AliESDtrack::kTRDout)!=0) && 
188         ((t->GetStatus()&AliESDtrack::kTRDStop)==0)){
189       AliTOFtrack *track = new AliTOFtrack(*t); 
190       track->SetSeedIndex(i);
191       t->UpdateTrackParams(track,AliESDtrack::kTOFout);    
192       new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
193       fNseedsTOF++;
194       delete track;
195     }
196
197     // Propagate the rest of TPCbp  
198
199     else {
200       AliTOFtrack *track = new AliTOFtrack(*t);
201       if(track->PropagateToInnerTOF(fHoles)){ // temporary solution
202         //      if(track->PropagateToInnerTOF(fGeom->GetHoles())){
203         track->SetSeedIndex(i);
204         t->UpdateTrackParams(track,AliESDtrack::kTOFout);    
205         new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
206         fNseedsTOF++;
207       }
208       delete track;
209     }
210   }
211
212   // Sort according uncertainties on track position 
213   fTracks->Sort();
214
215 }
216 //_________________________________________________________________________
217 void AliTOFtracker::MatchTracks( Bool_t mLastStep){
218
219   //Match ESD tracks to clusters in TOF
220
221   static const Double_t kMasses[]={
222     0.000511, 0.105658, 0.139570, 0.493677, 0.938272, 1.875613
223   };
224   
225   Int_t nSteps=(Int_t)(fTOFHeigth/0.1);
226
227   //PH Arrays (moved outside of the loop)
228   Float_t * trackPos[4];
229   for (Int_t ii=0; ii<4; ii++) trackPos[ii] = new Float_t[nSteps];
230   Int_t * clind[6];
231   for (Int_t ii=0;ii<6;ii++) clind[ii] = new Int_t[fN];
232   
233   for (Int_t i=0; i<fNseedsTOF; i++) {
234
235     AliTOFtrack *track =(AliTOFtrack*)fTracks->UncheckedAt(i);
236     AliESDtrack *t =(AliESDtrack*)fSeeds->UncheckedAt(track->GetSeedIndex());
237     if(t->GetTOFsignal()>0. ) continue;
238     AliTOFtrack *trackTOFin =new AliTOFtrack(*track);
239
240     // Some init 
241
242     Int_t         index[10000];
243     Float_t        dist[10000];
244     Float_t       cxpos[10000];
245     Float_t       crecL[10000];
246      
247     // Determine a window around the track
248
249     Double_t x,par[5]; 
250     trackTOFin->GetExternalParameters(x,par);
251     Double_t cov[15]; 
252     trackTOFin->GetExternalCovariance(cov);
253     Float_t scalefact=3.;    
254     Double_t dphi=
255       scalefact*
256       ((5*TMath::Sqrt(cov[0]) + 0.5*fDy + 2.5*TMath::Abs(par[2]))/fR); 
257     Double_t dz=
258       scalefact*
259       (5*TMath::Sqrt(cov[2]) + 0.5*fDz + 2.5*TMath::Abs(par[3]));
260     
261     Double_t phi=TMath::ATan2(par[0],x) + trackTOFin->GetAlpha();
262     if (phi<-TMath::Pi())phi+=2*TMath::Pi();
263     if (phi>=TMath::Pi())phi-=2*TMath::Pi();
264     Double_t z=par[1];   
265
266     Int_t nc=0;
267     
268     // find the clusters in the window of the track
269
270     for (Int_t k=FindClusterIndex(z-dz); k<fN; k++) {
271       AliTOFcluster *c=fClusters[k];
272       if (c->GetZ() > z+dz) break;
273       if (c->IsUsed()) continue;
274       
275       Double_t dph=TMath::Abs(c->GetPhi()-phi);
276       if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
277       if (TMath::Abs(dph)>dphi) continue;
278     
279       clind[0][nc] = c->GetDetInd(0);
280       clind[1][nc] = c->GetDetInd(1);
281       clind[2][nc] = c->GetDetInd(2);
282       clind[3][nc] = c->GetDetInd(3);
283       clind[4][nc] = c->GetDetInd(4);
284       clind[5][nc] = k;      
285       nc++;
286     }
287
288     //start fine propagation 
289
290     Int_t nStepsDone = 0;
291     for( Int_t istep=0; istep<nSteps; istep++){ 
292
293       Float_t xs=AliTOFGeometry::RinTOF()+istep*0.1;
294       Double_t ymax=xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
295
296       Bool_t skip=kFALSE;
297       Double_t ysect=trackTOFin->GetYat(xs,skip);
298       if(skip)break;
299       if (ysect > ymax) {
300         if (!trackTOFin->Rotate(AliTOFGeometry::GetAlpha())) {
301           break;
302         }
303       } else if (ysect <-ymax) {
304         if (!trackTOFin->Rotate(-AliTOFGeometry::GetAlpha())) {
305           break;
306         }
307       }
308
309       if(!trackTOFin->PropagateTo(xs)) {
310         break;
311       }
312
313       nStepsDone++;
314
315       // store the running point (Globalrf) - fine propagation     
316
317       Double_t x,y,z;
318       trackTOFin->GetGlobalXYZ(x,y,z);
319       trackPos[0][istep]= (Float_t) x;
320       trackPos[1][istep]= (Float_t) y;
321       trackPos[2][istep]= (Float_t) z;   
322       trackPos[3][istep]= trackTOFin->GetIntegratedLength();
323     }
324
325
326     Int_t nfound = 0;
327     for (Int_t istep=0; istep<nStepsDone; istep++) {
328
329       Bool_t isInside =kFALSE;
330       Float_t ctrackPos[3];     
331
332       ctrackPos[0]= trackPos[0][istep];
333       ctrackPos[1]= trackPos[1][istep];
334       ctrackPos[2]= trackPos[2][istep];
335
336       //now see whether the track matches any of the TOF clusters            
337
338       for (Int_t i=0; i<nc; i++){
339         Int_t cind[5];
340         cind[0]= clind[0][i];
341         cind[1]= clind[1][i];
342         cind[2]= clind[2][i];
343         cind[3]= clind[3][i];
344         cind[4]= clind[4][i];
345         Bool_t accept = kFALSE;
346         if( mLastStep)accept = (fGeom->DistanceToPad(cind,ctrackPos)<fdCut);
347         if(!mLastStep)accept = (fGeom->IsInsideThePad(cind,ctrackPos));
348         if(accept){
349           if(!mLastStep)isInside=kTRUE;
350           dist[nfound]=fGeom->DistanceToPad(cind,ctrackPos);
351           crecL[nfound]=trackPos[3][istep];
352           index[nfound]=clind[5][i]; // store cluster id            
353           cxpos[nfound]=AliTOFGeometry::RinTOF()+istep*0.1; //store prop.radius
354           nfound++;
355           if(isInside)break;
356         }//end if accept
357       } //end for on the clusters
358
359
360       if(isInside)break;
361     } //end for on the steps     
362
363
364
365     if (nfound == 0 ) {
366       fnunmatch++;
367       delete trackTOFin;
368       continue;
369     }
370     
371     fnmatch++;
372
373     // now choose the cluster to be matched with the track.
374
375     Int_t idclus=0;
376     Float_t  recL = 0.;
377     Float_t  xpos=0.;
378     Float_t  mindist=1000.;
379     for (Int_t iclus= 0; iclus<nfound;iclus++){
380       if (dist[iclus]< mindist){
381         mindist = dist[iclus];
382         xpos = cxpos[iclus];
383         idclus =index[iclus]; 
384         recL=crecL[iclus]+fDx*0.5;
385       }
386     }
387
388     AliTOFcluster *c=fClusters[idclus];
389     c->Use();
390
391     // Track length correction for matching Step 2 
392
393     if(mLastStep){
394       Float_t rc=TMath::Sqrt(c->GetR()*c->GetR() + c->GetZ()*c->GetZ());
395       Float_t rt=TMath::Sqrt(trackPos[0][70]*trackPos[0][70]
396                              +trackPos[1][70]*trackPos[1][70]
397                              +trackPos[2][70]*trackPos[2][70]);
398       Float_t dlt=rc-rt;      
399       recL=trackPos[3][70]+dlt;
400     }    
401
402     if (
403         (c->GetLabel(0)==TMath::Abs(trackTOFin->GetLabel()))
404         ||
405         (c->GetLabel(1)==TMath::Abs(trackTOFin->GetLabel()))
406         ||
407         (c->GetLabel(2)==TMath::Abs(trackTOFin->GetLabel()))
408         ) {
409       fngoodmatch++;
410     }
411     else{
412       fnbadmatch++;
413     }
414
415     delete trackTOFin;
416
417     Double_t tof=50*c->GetTDC()+32; // in ps
418     t->SetTOFsignal(tof);
419     t->SetTOFcluster(c->GetIndex());
420     Double_t time[10]; t->GetIntegratedTimes(time);
421     Double_t mom=t->GetP();
422     for(Int_t j=0;j<=5;j++){
423       Double_t mass=kMasses[j];
424       time[j]+=(recL-trackPos[3][0])/3e-2*TMath::Sqrt(mom*mom+mass*mass)/mom;
425     }
426
427     AliTOFtrack *trackTOFout = new AliTOFtrack(*t); 
428     trackTOFout->PropagateTo(xpos);
429     t->UpdateTrackParams(trackTOFout,AliESDtrack::kTOFout);    
430     t->SetIntegratedLength(recL);
431     t->SetIntegratedTimes(time);
432
433     delete trackTOFout;
434   }
435   for (Int_t ii=0; ii<4; ii++) delete [] trackPos[ii];
436   for (Int_t ii=0;ii<6;ii++) delete [] clind[ii];
437 }
438 //_________________________________________________________________________
439 Int_t AliTOFtracker::LoadClusters(TTree *dTree) {
440   //--------------------------------------------------------------------
441   //This function loads the TOF clusters
442   //--------------------------------------------------------------------
443
444   TBranch *branch=dTree->GetBranch("TOF");
445   if (!branch) { 
446     Error("LoadClusters"," can't get the branch with the TOF digits !\n");
447     return 1;
448   }
449
450   TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy;
451   branch->SetAddress(&digits);
452
453   dTree->GetEvent(0);
454   Int_t nd=digits->GetEntriesFast();
455   Info("LoadClusters","number of digits: %d",nd);
456
457   for (Int_t i=0; i<nd; i++) {
458     AliTOFdigit *d=(AliTOFdigit*)digits->UncheckedAt(i);
459     Int_t dig[5]; Float_t g[3];
460     dig[0]=d->GetSector();
461     dig[1]=d->GetPlate();
462     dig[2]=d->GetStrip();
463     dig[3]=d->GetPadz();
464     dig[4]=d->GetPadx();
465
466     fGeom->GetPos(dig,g);
467
468     Double_t h[5];
469     h[0]=TMath::Sqrt(g[0]*g[0]+g[1]*g[1]);
470     h[1]=TMath::ATan2(g[1],g[0]); h[2]=g[2]; 
471     h[3]=d->GetTdc(); h[4]=d->GetAdc();
472
473     AliTOFcluster *cl=new AliTOFcluster(h,d->GetTracks(),dig,i);
474     InsertCluster(cl);
475   }  
476
477   return 0;
478 }
479 //_________________________________________________________________________
480 void AliTOFtracker::UnloadClusters() {
481   //--------------------------------------------------------------------
482   //This function unloads TOF clusters
483   //--------------------------------------------------------------------
484   for (Int_t i=0; i<fN; i++) delete fClusters[i];
485   fN=0;
486 }
487
488 //_________________________________________________________________________
489 Int_t AliTOFtracker::InsertCluster(AliTOFcluster *c) {
490   //--------------------------------------------------------------------
491   //This function adds a cluster to the array of clusters sorted in Z
492   //--------------------------------------------------------------------
493   if (fN==kMaxCluster) {
494     Error("InsertCluster","Too many clusters !\n");
495     return 1;
496   }
497
498   if (fN==0) {fClusters[fN++]=c; return 0;}
499   Int_t i=FindClusterIndex(c->GetZ());
500   memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTOFcluster*));
501   fClusters[i]=c; fN++;
502
503   return 0;
504 }
505
506 //_________________________________________________________________________
507 Int_t AliTOFtracker::FindClusterIndex(Double_t z) const {
508   //--------------------------------------------------------------------
509   // This function returns the index of the nearest cluster 
510   //--------------------------------------------------------------------
511   if (fN==0) return 0;
512   if (z <= fClusters[0]->GetZ()) return 0;
513   if (z > fClusters[fN-1]->GetZ()) return fN;
514   Int_t b=0, e=fN-1, m=(b+e)/2;
515   for (; b<e; m=(b+e)/2) {
516     if (z > fClusters[m]->GetZ()) b=m+1;
517     else e=m; 
518   }
519   return m;
520 }
521