uuhhhhaaa - what did I change? Added the run number and fixed coding conventions ;-)
[u/mrichter/AliRoot.git] / TOF / AliTOFtrackerMI.cxx
CommitLineData
d88fbf15 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 **************************************************************************/
0e46b9ae 15
16//-----------------------------------------------------------------//
17// //
18// AliTOFtracker Class //
19// Task: Perform association of the ESD tracks to TOF Clusters //
20// and Update ESD track with associated TOF Cluster parameters //
21// //
22//-----------------------------------------------------------------//
d88fbf15 23
24#include <Rtypes.h>
0e46b9ae 25
d88fbf15 26#include "TClonesArray.h"
0e46b9ae 27#include "TGeoManager.h"
28#include "TTree.h"
d88fbf15 29#include "TTreeStream.h"
0e46b9ae 30
31#include "AliRun.h"
32#include "AliESD.h"
33#include "AliESDtrack.h"
34
a533f541 35#include "AliTOFcalib.h"
0e46b9ae 36#include "AliTOFcluster.h"
37#include "AliTOFGeometry.h"
38#include "AliTOFtrackerMI.h"
39#include "AliTOFtrack.h"
40
41extern TGeoManager *gGeoManager;
d88fbf15 42
43ClassImp(AliTOFtrackerMI)
44
45//_____________________________________________________________________________
58d8d9a3 46AliTOFtrackerMI::AliTOFtrackerMI(AliTOFGeometry * geom, Double_t parPID[2]):
47 fGeom(geom),
48 fTOFpid(new AliTOFpidESD(parPID)),
49 fHoles(kFALSE),
50 fN(0),
51 fNseeds(0),
52 fNseedsTOF(0),
53 fngoodmatch(0),
54 fnbadmatch(0),
55 fnunmatch(0),
56 fnmatch(0),
57 fR(378.),
58 fTOFHeigth(15.3),
59 fdCut(3.),
60 fDx(1.5),
61 fDy(0),
62 fDz(0),
63 fTracks(0x0),
64 fSeeds(0x0),
65 fDebugStreamer(0x0)
66 {
d88fbf15 67 //AliTOFtrackerMI main Ctor
68
d3c7bfac 69 //fHoles=true;
d88fbf15 70 fDy=AliTOFGeometry::XPad();
71 fDz=AliTOFGeometry::ZPad();
d88fbf15 72 fDebugStreamer = new TTreeSRedirector("TOFdebug.root");
d3c7bfac 73 //Init(); // temporary solution to know about Holes/no Holes
74 fHoles=geom->GetHoles();
d88fbf15 75}
76//_____________________________________________________________________________
58d8d9a3 77AliTOFtrackerMI::AliTOFtrackerMI(const AliTOFtrackerMI &t):
78 AliTracker(),
79 fGeom(0x0),
80 fTOFpid(0x0),
81 fHoles(kFALSE),
82 fN(0),
83 fNseeds(0),
84 fNseedsTOF(0),
85 fngoodmatch(0),
86 fnbadmatch(0),
87 fnunmatch(0),
88 fnmatch(0),
89 fR(378.),
90 fTOFHeigth(15.3),
91 fdCut(3.),
92 fDx(1.5),
93 fDy(0),
94 fDz(0),
95 fTracks(0x0),
96 fSeeds(0x0),
97 fDebugStreamer(0x0)
98 {
d88fbf15 99 //AliTOFtrackerMI copy Ctor
100
101 fHoles=t.fHoles;
102 fNseeds=t.fNseeds;
103 fNseedsTOF=t.fNseedsTOF;
104 fngoodmatch=t.fngoodmatch;
105 fnbadmatch=t.fnbadmatch;
106 fnunmatch=t.fnunmatch;
107 fnmatch=t.fnmatch;
108 fGeom = t.fGeom;
109 fTOFpid = t.fTOFpid;
110 fR=t.fR;
111 fTOFHeigth=t.fTOFHeigth;
112 fdCut=t.fdCut;
113 fDy=t.fDy;
114 fDz=t.fDz;
115 fDx=1.5;
116 fSeeds=t.fSeeds;
117 fTracks=t.fTracks;
118 fN=t.fN;
119}
120
121//_________________________________________________________________________________
7aeeaf38 122AliTOFtrackerMI& AliTOFtrackerMI::operator=(const AliTOFtrackerMI &t)
123{
124 //AliTOFtrackerMI assignment operator
125
126 this->fHoles=t.fHoles;
127 this->fNseeds=t.fNseeds;
128 this->fNseedsTOF=t.fNseedsTOF;
129 this->fngoodmatch=t.fngoodmatch;
130 this->fnbadmatch=t.fnbadmatch;
131 this->fnunmatch=t.fnunmatch;
132 this->fnmatch=t.fnmatch;
133 this->fGeom = t.fGeom;
134 this->fTOFpid = t.fTOFpid;
135 this->fR=t.fR;
136 this->fTOFHeigth=t.fTOFHeigth;
137 this->fdCut=t.fdCut;
138 this->fDy=t.fDy;
139 this->fDz=t.fDz;
140 this->fDx=t.fDx;
141 this->fSeeds=t.fSeeds;
142 this->fTracks=t.fTracks;
143 this->fN=t.fN;
144 return *this;
145
146}
147
148//_____________________________________________________________________________
d88fbf15 149AliTOFtrackerMI::~AliTOFtrackerMI(){
150 //
151 //
152 //
153 if (fDebugStreamer) {
154 //fDebugStreamer->Close();
155 delete fDebugStreamer;
156 }
157 delete fTOFpid;
158}
159
160//_____________________________________________________________________________
d3c7bfac 161/*
d88fbf15 162void AliTOFtrackerMI::Init() {
163
164// temporary solution to know about Holes/no Holes, will be implemented as
165// an AliTOFGeometry getter
166
167 AliModule* frame=gAlice->GetModule("FRAME");
168
169 if(!frame) {
0e46b9ae 170 AliError("Could Not load FRAME! Assume Frame with Holes ");
d88fbf15 171 fHoles=true;
172 } else{
173 if(frame->IsVersion()==1) {fHoles=false;}
174 else {fHoles=true;}
175 }
176}
d3c7bfac 177*/
d88fbf15 178//_____________________________________________________________________________
179Int_t AliTOFtrackerMI::PropagateBack(AliESD* event) {
180 //
181 // Gets seeds from ESD event and Match with TOF Clusters
182 //
183
184
185 //Initialise some counters
186
187 fNseeds=0;
188 fNseedsTOF=0;
189 fngoodmatch=0;
190 fnbadmatch=0;
191 fnunmatch=0;
192 fnmatch=0;
193
194 Int_t ntrk=event->GetNumberOfTracks();
195 fNseeds = ntrk;
196 fSeeds= new TClonesArray("AliESDtrack");
197 TClonesArray &aESDTrack = *fSeeds;
198
199
200 //Load ESD tracks into a local Array of ESD Seeds
201
202 for (Int_t i=0; i<fNseeds; i++) {
203 AliESDtrack *t=event->GetTrack(i);
204 new(aESDTrack[i]) AliESDtrack(*t);
205 }
206
207 //Prepare ESD tracks candidates for TOF Matching
208 CollectESD();
209
210 //First Step with Strict Matching Criterion
211 //MatchTracks(kFALSE);
212
213 //Second Step with Looser Matching Criterion
214 //MatchTracks(kTRUE);
215 MatchTracksMI(kFALSE); // assign track to clusters
216 MatchTracksMI(kTRUE); // assign clusters to esd
217
218 Info("PropagateBack","Number of matched tracks: %d",fnmatch);
219 Info("PropagateBack","Number of good matched tracks: %d",fngoodmatch);
220 Info("PropagateBack","Number of bad matched tracks: %d",fnbadmatch);
221
222 //Update the matched ESD tracks
223
224 for (Int_t i=0; i<ntrk; i++) {
225 AliESDtrack *t=event->GetTrack(i);
226 AliESDtrack *seed =(AliESDtrack*)fSeeds->UncheckedAt(i);
227 if(seed->GetTOFsignal()>0){
228 t->SetTOFsignal(seed->GetTOFsignal());
229 t->SetTOFcluster(seed->GetTOFcluster());
230 Int_t tlab[3];
231 seed->GetTOFLabel(tlab);
232 t->SetTOFLabel(tlab);
233 AliTOFtrack *track = new AliTOFtrack(*seed);
234 Float_t info[10];
235 Double_t times[10];
236 seed->GetTOFInfo(info);
237 seed->GetIntegratedTimes(times);
238 t->UpdateTrackParams(track,AliESDtrack::kTOFout);
239 t->SetIntegratedLength(seed->GetIntegratedLength());
240 t->SetIntegratedTimes(times);
a533f541 241 t->SetTOFsignalToT(seed->GetTOFsignalToT());
242 t->SetTOFCalChannel(seed->GetTOFCalChannel());
d88fbf15 243 //
244 t->SetTOFInfo(info);
245 delete track;
246 }
247 }
248
249
250 //Make TOF PID
251 fTOFpid->MakePID(event);
252
253 if (fSeeds) {
254 fSeeds->Delete();
255 delete fSeeds;
256 fSeeds = 0x0;
257 }
258 if (fTracks) {
259 fTracks->Delete();
260 delete fTracks;
261 fTracks = 0x0;
262 }
263 return 0;
264
265}
266//_________________________________________________________________________
267void AliTOFtrackerMI::CollectESD() {
268 //prepare the set of ESD tracks to be matched to clusters in TOF
269
270 fTracks= new TClonesArray("AliTOFtrack");
271 TClonesArray &aTOFTrack = *fTracks;
272 Int_t c0=0;
273 Int_t c1=0;
274 for (Int_t i=0; i<fNseeds; i++) {
275
276 AliESDtrack *t =(AliESDtrack*)fSeeds->UncheckedAt(i);
277 if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
278
279 // TRD good tracks, already propagated at 371 cm
280
281 AliTOFtrack *track = new AliTOFtrack(*t); // New
282 Double_t x = track->GetX(); //New
283
284 if (((t->GetStatus()&AliESDtrack::kTRDout)!=0 ) &&
d3c7bfac 285 ( x >= fGeom->RinTOF()) ){
d88fbf15 286 track->SetSeedIndex(i);
287 t->UpdateTrackParams(track,AliESDtrack::kTOFout);
288 new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
289 fNseedsTOF++;
290 c0++;
291 delete track;
292 }
293
294 // Propagate the rest of TPCbp
295
296 else {
297 if(track->PropagateToInnerTOF(fHoles)){ // temporary solution
298 // if(track->PropagateToInnerTOF(fGeom->GetHoles())){
299 track->SetSeedIndex(i);
300 t->UpdateTrackParams(track,AliESDtrack::kTOFout);
301 new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
302 fNseedsTOF++;
303 c1++;
304 }
305 delete track;
306 }
307 }
308 //
309 //
310 printf("TRD\tOn\t%d\tOff\t%d\n",c0,c1);
311
312 // Sort according uncertainties on track position
313 fTracks->Sort();
314
315}
316
317
318
319
320
321
322
323//
324//
325//_________________________________________________________________________
326void AliTOFtrackerMI::MatchTracks( Bool_t /*mLastStep*/){
327 return;
328}
329//
330//
331//_________________________________________________________________________
332void AliTOFtrackerMI::MatchTracksMI(Bool_t mLastStep){
333
334 //Match ESD tracks to clusters in TOF
335 const Float_t kTofOffset = 26; // time offset
336 const Float_t kMinQuality = -6.; // minimal quality
337 const Float_t kMaxQualityD = 1.; // max delta quality if cluster used
338 const Float_t kForbiddenR = 0.1; // minimal PID according TPC
339
340 static const Double_t kMasses[]={
341 0.000511, 0.105658, 0.139570, 0.493677, 0.938272, 1.875613
342 };
343
344 Int_t nSteps=(Int_t)(fTOFHeigth/0.1);
345
a533f541 346 AliTOFcalib *calib = new AliTOFcalib(fGeom);
347
d88fbf15 348 //PH Arrays (moved outside of the loop)
349 Float_t * trackPos[4];
350 for (Int_t ii=0; ii<4; ii++) trackPos[ii] = new Float_t[nSteps];
351 Int_t * clind[6];
352 for (Int_t ii=0;ii<6;ii++) clind[ii] = new Int_t[fN];
353
354 // Some init
355
356 Int_t index[1000];
357 Float_t dist3D[1000][6];
358 Double_t times[1000][6];
359 Float_t mintimedist[1000];
360 Float_t likelihood[1000];
361 Float_t length[1000];
362 AliTOFcluster *clusters[1000];
363 Double_t tpcpid[5];
364 dist3D[0][0]=1;
365
366 for (Int_t i=0; i<fNseedsTOF; i++) {
367
368 AliTOFtrack *track =(AliTOFtrack*)fTracks->UncheckedAt(i);
369 AliESDtrack *t =(AliESDtrack*)fSeeds->UncheckedAt(track->GetSeedIndex());
370 Bool_t hasTime = ( (t->GetStatus()& AliESDtrack::kTIME)>0) ? kTRUE:kFALSE; // did we integrate time
371 Float_t trdquality = t->GetTRDQuality();
372 //
373 // Normalize tpc pid
374 //
375 t->GetTPCpid(tpcpid);
376 Double_t sumpid=0;
377 for (Int_t ipid=0;ipid<5;ipid++){
378 sumpid+=tpcpid[ipid];
379 }
380 for (Int_t ipid=0;ipid<5;ipid++){
381 if (sumpid>0) tpcpid[ipid]/=sumpid;
382 else{
383 tpcpid[ipid]=0.2;
384 }
385 }
386
387 if (trdquality<0) continue; // no chance
388 //
389 AliTOFtrack *trackTOFin =new AliTOFtrack(*track);
390 //
391 //propagat track to the middle of TOF
392 //
d3c7bfac 393 Float_t xs = 378.2; // should be defined in the TOF geometry
d88fbf15 394 Double_t ymax=xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
d3c7bfac 395 Bool_t skip=kFALSE;
d88fbf15 396 Double_t ysect=trackTOFin->GetYat(xs,skip);
397 if (skip){
398 xs = 372.;
399 ymax=xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
400 ysect=trackTOFin->GetYat(xs,skip);
401 }
402 if (ysect > ymax) {
403 if (!trackTOFin->Rotate(AliTOFGeometry::GetAlpha())) {
404 continue;
405 }
406 } else if (ysect <-ymax) {
407 if (!trackTOFin->Rotate(-AliTOFGeometry::GetAlpha())) {
408 continue;
409 }
410 }
411 if(!trackTOFin->PropagateTo(xs)) {
412 continue;
413 }
414 //
415 // Determine a window around the track
416 //
417 Double_t x,par[5];
418 trackTOFin->GetExternalParameters(x,par);
419 Double_t cov[15];
420 trackTOFin->GetExternalCovariance(cov);
421 Float_t scalefact=3.;
422 Double_t dphi=
423 scalefact*
424 ((5*TMath::Sqrt(cov[0]) + 3.*fDy +10.*TMath::Abs(par[2]))/fR);
425 Double_t dz=
426 scalefact*
427 (5*TMath::Sqrt(cov[2]) + 3.*fDz +10.*TMath::Abs(par[3]));
428
429 Double_t phi=TMath::ATan2(par[0],x) + trackTOFin->GetAlpha();
430 if (phi<-TMath::Pi())phi+=2*TMath::Pi();
431 if (phi>=TMath::Pi())phi-=2*TMath::Pi();
432 Double_t z=par[1];
433
434 Int_t nc =0;
435 Int_t nfound =0;
436 // find the clusters in the window of the track
437
438 for (Int_t k=FindClusterIndex(z-dz); k<fN; k++) {
439 AliTOFcluster *c=fClusters[k];
440 if (c->GetZ() > z+dz) break;
441 // if (c->IsUsed()) continue;
442
443 Double_t dph=TMath::Abs(c->GetPhi()-phi);
444 if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
445 if (TMath::Abs(dph)>dphi) continue;
446
447 clind[0][nc] = c->GetDetInd(0);
448 clind[1][nc] = c->GetDetInd(1);
449 clind[2][nc] = c->GetDetInd(2);
450 clind[3][nc] = c->GetDetInd(3);
451 clind[4][nc] = c->GetDetInd(4);
452 clind[5][nc] = k;
453 nc++;
454 }
455
456 //
457 // select close clusters
458 //
459 Double_t mom=t->GetP();
460 // Bool_t dump = kTRUE;
461 for (Int_t icl=0; icl<nc; icl++){
462 Float_t distances[5];
463 if (nfound>=1000) break;
464 index[nfound]=clind[5][icl];
465 AliTOFcluster *cluster = fClusters[clind[5][icl]];
466 GetLinearDistances(trackTOFin,cluster, distances);
467 dist3D[nfound][0] = distances[4];
468 dist3D[nfound][1] = distances[1];
469 dist3D[nfound][2] = distances[2];
470 // cut on distance
471 if (TMath::Abs(dist3D[nfound][1])>20 || TMath::Abs(dist3D[nfound][2])>20) continue;
472 //
473 GetLikelihood(distances[1],distances[2],cov,trackTOFin, dist3D[nfound][3], dist3D[nfound][4]);
474 //
475 // cut on likelihood
476 if (dist3D[nfound][3]*dist3D[nfound][4]<0.00000000000001) continue; // log likelihood
477 if (TMath::Log(dist3D[nfound][3]*dist3D[nfound][4])<-9) continue; // log likelihood
478 //
479 clusters[nfound] = cluster;
480 //
481 //length and TOF updates
482 trackTOFin->GetIntegratedTimes(times[nfound]);
483 length[nfound] = trackTOFin->GetIntegratedLength();
484 length[nfound]+=distances[4];
485 mintimedist[nfound]=1000;
486 Double_t tof2=AliTOFGeometry::TdcBinWidth()*cluster->GetTDC()+kTofOffset; // in ps
487 // Float_t tgamma = TMath::Sqrt(cluster->GetR()*cluster->GetR()+cluster->GetZ()*cluster->GetZ())/0.03; //time for "primary" gamma
488 //if (trackTOFin->GetPt()<0.7 && TMath::Abs(tgamma-tof2)<350) continue; // gamma conversion candidate - TEMPORARY
489 for(Int_t j=0;j<=5;j++){
490 Double_t mass=kMasses[j];
491 times[nfound][j]+=distances[4]/3e-2*TMath::Sqrt(mom*mom+mass*mass)/mom; // add time distance
492 if ( TMath::Abs(times[nfound][j]-tof2)<mintimedist[nfound] && tpcpid[j]>kForbiddenR){
493 mintimedist[nfound]=TMath::Abs(times[nfound][j]-tof2);
494 }
495 }
496 //
497 Float_t liketime = TMath::Exp(-mintimedist[nfound]/90.)+0.25*TMath::Exp(-mintimedist[nfound]/180.);
498 if (!hasTime) liketime=0.2;
499 likelihood[nfound] = TMath::Log(dist3D[nfound][3]*dist3D[nfound][4]*liketime);
500
501 if (TMath::Log(dist3D[nfound][3]*dist3D[nfound][4])<-1){
502 if (likelihood[nfound]<-9.) continue;
503 }
504 //
505 nfound++;
506 }
507 if (nfound == 0 ) {
508 fnunmatch++;
509 delete trackTOFin;
510 continue;
511 }
512 //
513 //choose the best cluster
514 //
515 Float_t quality[1000];
516 Int_t index[1000];
517 //
518 AliTOFcluster * cgold=0;
519 Int_t igold =-1;
520 for (Int_t icl=0;icl<nfound;icl++){
521 quality[icl] = dist3D[icl][3]*dist3D[icl][4];
522 }
523 TMath::Sort(nfound,quality,index,kTRUE);
524 igold = index[0];
525 cgold = clusters[igold];
526 if (nfound>1 &&likelihood[index[1]]>likelihood[index[0]]){
527 if ( quality[index[0]]<quality[index[1]]+0.5){
528 igold = index[1];
529 cgold = clusters[igold];
530 }
531 }
532 //
533 //
534 Float_t qualityGold = TMath::Log(0.0000001+(quality[igold]*(0.1+TMath::Exp(-mintimedist[igold]/80.))*(0.2+trdquality)));
535 if (!mLastStep){
536 if (cgold->GetQuality()<qualityGold) cgold->SetQuality(qualityGold);
537 continue;
538 }
539 //
540 if (mLastStep){
541 //signed better cluster
542 if (cgold->GetQuality()>qualityGold+kMaxQualityD) continue;
543 if (2.*qualityGold-cgold->GetQuality()<kMinQuality) continue;
544 }
545
546 Int_t inonfake=-1;
547
548 for (Int_t icl=0;icl<nfound;icl++){
549 if (
550 (clusters[index[icl]]->GetLabel(0)==TMath::Abs(trackTOFin->GetLabel()))
551 ||
552 (clusters[index[icl]]->GetLabel(1)==TMath::Abs(trackTOFin->GetLabel()))
553 ||
554 (clusters[index[icl]]->GetLabel(2)==TMath::Abs(trackTOFin->GetLabel()))
555 ) {
556 inonfake = icl;
557 break;
558 }
559 }
560 fnmatch++;;
561 if (inonfake==0) fngoodmatch++;
562 else{
563 fnbadmatch++;
564 }
565
566 Int_t tlab[3];
567 tlab[0]=cgold->GetLabel(0);
568 tlab[1]=cgold->GetLabel(1);
569 tlab[2]=cgold->GetLabel(2);
570 // Double_t tof2=25.*cgold->GetTDC()-350; // in ps
571 Double_t tof2=AliTOFGeometry::TdcBinWidth()*cgold->GetTDC()+kTofOffset; // in ps
572 Float_t tgamma = TMath::Sqrt(cgold->GetR()*cgold->GetR()+cgold->GetZ()*cgold->GetZ())/0.03;
573 Float_t info[11]={dist3D[igold][0],dist3D[igold][1],dist3D[igold][2],dist3D[igold][3],dist3D[igold][4],mintimedist[igold],
574 -1,tgamma, qualityGold,cgold->GetQuality(),0};
575 // GetLinearDistances(trackTOFin,cgold,&info[6]);
576 if (inonfake>=0){
577 info[6] = inonfake;
578 // info[7] = mintimedist[index[inonfake]];
579 }
580 //
a533f541 581 // Store quantities to be used for TOF Calibration
7aeeaf38 582 Float_t tToT=cgold->GetToT(); // in ps
583 t->SetTOFsignalToT(tToT);
a533f541 584 Int_t ind[5];
585 ind[0]=cgold->GetDetInd(0);
586 ind[1]=cgold->GetDetInd(1);
587 ind[2]=cgold->GetDetInd(2);
588 ind[3]=cgold->GetDetInd(3);
589 ind[4]=cgold->GetDetInd(4);
590 Int_t calindex = calib->GetIndex(ind);
591 t->SetTOFCalChannel(calindex);
592
d88fbf15 593 t->SetTOFInfo(info);
594 t->SetTOFsignal(tof2);
595 t->SetTOFcluster(cgold->GetIndex());
596 AliTOFtrack *trackTOFout = new AliTOFtrack(*t);
597 trackTOFout->PropagateTo(378.);
598 t->UpdateTrackParams(trackTOFout,AliESDtrack::kTOFout);
599 t->SetIntegratedLength(length[igold]);
600 t->SetIntegratedTimes(times[igold]);
601 t->SetTOFLabel(tlab);
602 //
603 delete trackTOFin;
604 delete trackTOFout;
605 //
606 }
607 //
608 //
609 //
610 for (Int_t ii=0; ii<4; ii++) delete [] trackPos[ii];
611 for (Int_t ii=0;ii<6;ii++) delete [] clind[ii];
a533f541 612 delete calib;
d88fbf15 613}
614
615
616
617
618// //_________________________________________________________________________
619// Int_t AliTOFtrackerMI::LoadClusters(TTree *dTree) {
620// //--------------------------------------------------------------------
621// //This function loads the TOF clusters
622// //--------------------------------------------------------------------
623
624// TBranch *branch=dTree->GetBranch("TOF");
625// if (!branch) {
0e46b9ae 626// AliError(" can't get the branch with the TOF digits !");
d88fbf15 627// return 1;
628// }
629
630// TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy;
631// branch->SetAddress(&digits);
632
633// dTree->GetEvent(0);
634// Int_t nd=digits->GetEntriesFast();
635// Info("LoadClusters","number of digits: %d",nd);
636
637// for (Int_t i=0; i<nd; i++) {
638// AliTOFdigit *d=(AliTOFdigit*)digits->UncheckedAt(i);
639// Int_t dig[5]; Float_t g[3];
640// dig[0]=d->GetSector();
641// dig[1]=d->GetPlate();
642// dig[2]=d->GetStrip();
643// dig[3]=d->GetPadz();
644// dig[4]=d->GetPadx();
645
646// fGeom->GetPos(dig,g);
647
648// Double_t h[5];
649// h[0]=TMath::Sqrt(g[0]*g[0]+g[1]*g[1]);
650// h[1]=TMath::ATan2(g[1],g[0]); h[2]=g[2];
651// h[3]=d->GetTdc(); h[4]=d->GetAdc();
652
653// AliTOFcluster *cl=new AliTOFcluster(h,d->GetTracks(),dig,i);
654// InsertCluster(cl);
655// }
656
657// return 0;
658// }
659// //_________________________________________________________________________
660// void AliTOFtrackerMI::UnloadClusters() {
661// //--------------------------------------------------------------------
662// //This function unloads TOF clusters
663// //--------------------------------------------------------------------
664// for (Int_t i=0; i<fN; i++) delete fClusters[i];
665// fN=0;
666// }
667
668
669
670//_________________________________________________________________________
671Int_t AliTOFtrackerMI::LoadClusters(TTree *cTree) {
672 //--------------------------------------------------------------------
673 //This function loads the TOF clusters
674 //--------------------------------------------------------------------
675
676 TBranch *branch=cTree->GetBranch("TOF");
677 if (!branch) {
678 AliError("can't get the branch with the TOF clusters !");
679 return 1;
680 }
681
682 TClonesArray dummy("AliTOFcluster",10000), *clusters=&dummy;
683 branch->SetAddress(&clusters);
684
685 cTree->GetEvent(0);
686 Int_t nc=clusters->GetEntriesFast();
687 AliInfo(Form("Number of clusters: %d",nc));
688
689 for (Int_t i=0; i<nc; i++) {
690 AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
691 /*
692 for (Int_t jj=0; jj<5; jj++) dig[jj]=c->GetDetInd(jj);
693
694 Double_t h[5];
695 h[0]=c->GetR();
696 h[1]=c->GetPhi();
697 h[2]=c->GetZ();
698 h[3]=c->GetTDC();
699 h[4]=c->GetADC();
700
701 Int_t indexDig[3];
702 for (Int_t jj=0; jj<3; jj++) indexDig[jj] = c->GetLabel(jj);
703
704 AliTOFcluster *cl=new AliTOFcluster(h,c->GetTracks(),dig,i);
705 */
706
707 // fClusters[i]=c; fN++;
708 fClusters[i]=new AliTOFcluster(*c); fN++;
709
710 //AliInfo(Form("%4i %4i %f %f %f %f %f %2i %1i %2i %1i %2i",i, fClusters[i]->GetIndex(),fClusters[i]->GetZ(),fClusters[i]->GetR(),fClusters[i]->GetPhi(), fClusters[i]->GetTDC(),fClusters[i]->GetADC(),fClusters[i]->GetDetInd(0),fClusters[i]->GetDetInd(1),fClusters[i]->GetDetInd(2),fClusters[i]->GetDetInd(3),fClusters[i]->GetDetInd(4)));
711 //AliInfo(Form("%i %f",i, fClusters[i]->GetZ()));
712 }
713
714 //AliInfo(Form("Number of clusters: %d",fN));
715
716 return 0;
717}
718//_________________________________________________________________________
719void AliTOFtrackerMI::UnloadClusters() {
720 //--------------------------------------------------------------------
721 //This function unloads TOF clusters
722 //--------------------------------------------------------------------
723 for (Int_t i=0; i<fN; i++) {
724 delete fClusters[i];
725 fClusters[i] = 0x0;
726 }
727 fN=0;
728}
729
730
731
732
733//_________________________________________________________________________
734Int_t AliTOFtrackerMI::InsertCluster(AliTOFcluster *c) {
735 //--------------------------------------------------------------------
736 //This function adds a cluster to the array of clusters sorted in Z
737 //--------------------------------------------------------------------
738 if (fN==kMaxCluster) {
0e46b9ae 739 AliError("Too many clusters !");
d88fbf15 740 return 1;
741 }
742
743 if (fN==0) {fClusters[fN++]=c; return 0;}
744 Int_t i=FindClusterIndex(c->GetZ());
745 memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTOFcluster*));
746 fClusters[i]=c; fN++;
747
748 return 0;
749}
750
751//_________________________________________________________________________
752Int_t AliTOFtrackerMI::FindClusterIndex(Double_t z) const {
753 //--------------------------------------------------------------------
754 // This function returns the index of the nearest cluster
755 //--------------------------------------------------------------------
756 if (fN==0) return 0;
757 if (z <= fClusters[0]->GetZ()) return 0;
758 if (z > fClusters[fN-1]->GetZ()) return fN;
759 Int_t b=0, e=fN-1, m=(b+e)/2;
760 for (; b<e; m=(b+e)/2) {
761 if (z > fClusters[m]->GetZ()) b=m+1;
762 else e=m;
763 }
764 return m;
765}
766
767
768
769Float_t AliTOFtrackerMI::GetLinearDistances(AliTOFtrack * track, AliTOFcluster *cluster, Float_t distances[5])
770{
771 //
772 // calclates distance between cluster and track
773 // use linear aproximation
774 //
775 const Float_t kRaddeg = 180/3.14159265358979312;
776 //
777 // Float_t tiltangle = fGeom->GetAngles(cluster->fdetIndex[1],cluster->fdetIndex[2])/kRaddeg; //tiltangle
778 Int_t cind[5];
779 cind[0]= cluster->GetDetInd(0);
780 cind[1]= cluster->GetDetInd(1);
781 cind[2]= cluster->GetDetInd(2);
782 cind[3]= cluster->GetDetInd(3);
783 cind[4]= cluster->GetDetInd(4);
784 Float_t tiltangle = fGeom->GetAngles(cluster->GetDetInd(1),cluster->GetDetInd(2))/kRaddeg; //tiltangle
785
786 Float_t cpos[3]; //cluster position
787 Float_t cpos0[3]; //cluster position
788 // fGeom->GetPos(cluster->fdetIndex,cpos);
789 //fGeom->GetPos(cluster->fdetIndex,cpos0);
790 //
791 fGeom->GetPos(cind,cpos);
792 fGeom->GetPos(cind,cpos0);
793
794 Float_t phi = TMath::ATan2(cpos[1],cpos[0]);
795 if(phi<0) phi=2.*TMath::Pi()+phi;
796 // Get the local angle in the sector philoc
797 Float_t phiangle = (Int_t (phi*kRaddeg/20.) + 0.5)*20./kRaddeg;
798 //
799 Double_t v0[3];
800 Double_t dir[3];
6c94f330 801 track->GetXYZ(v0);
802 track->GetPxPyPz(dir);
d88fbf15 803 dir[0]/=track->GetP();
804 dir[1]/=track->GetP();
805 dir[2]/=track->GetP();
806 //
807 //
808 //rotate 0
809 Float_t sinphi = TMath::Sin(phiangle);
810 Float_t cosphi = TMath::Cos(phiangle);
811 Float_t sinth = TMath::Sin(tiltangle);
812 Float_t costh = TMath::Cos(tiltangle);
813 //
814 Float_t temp;
815 temp = cpos[0]*cosphi+cpos[1]*sinphi;
816 cpos[1] = -cpos[0]*sinphi+cpos[1]*cosphi;
817 cpos[0] = temp;
818 temp = v0[0]*cosphi+v0[1]*sinphi;
819 v0[1] = -v0[0]*sinphi+v0[1]*cosphi;
820 v0[0] = temp;
821 //
822 temp = cpos[0]*costh+cpos[2]*sinth;
823 cpos[2] = -cpos[0]*sinth+cpos[2]*costh;
824 cpos[0] = temp;
825 temp = v0[0]*costh+v0[2]*sinth;
826 v0[2] = -v0[0]*sinth+v0[2]*costh;
827 v0[0] = temp;
828 //
829 //
830 //rotate direction vector
831 //
832 temp = dir[0]*cosphi+dir[1]*sinphi;
833 dir[1] = -dir[0]*sinphi+dir[1]*cosphi;
834 dir[0] = temp;
835 //
836 temp = dir[0]*costh+dir[2]*sinth;
837 dir[2] = -dir[0]*sinth+dir[2]*costh;
838 dir[0] = temp;
839 //
840 Float_t v3[3];
841 Float_t k = (cpos[0]-v0[0])/dir[0];
842 v3[0] = v0[0]+k*dir[0];
843 v3[1] = v0[1]+k*dir[1];
844 v3[2] = v0[2]+k*dir[2];
845 //
846 distances[0] = v3[0]-cpos[0];
847 distances[1] = v3[1]-cpos[1];
848 distances[2] = v3[2]-cpos[2];
849 distances[3] = TMath::Sqrt( distances[0]*distances[0]+distances[1]*distances[1]+distances[2]*distances[2]); //distance
850 distances[4] = k; //length
851
852 //
853 // Debuging part of the matching
854 //
855 if (track->GetLabel()==cluster->GetLabel(0) ||track->GetLabel()==cluster->GetLabel(1) ){
856 TTreeSRedirector& cstream = *fDebugStreamer;
857 Float_t tdc = cluster->GetTDC();
858 cstream<<"Tracks"<<
859 "TOF.="<<track<<
860 "Cx="<<cpos0[0]<<
861 "Cy="<<cpos0[1]<<
862 "Cz="<<cpos0[2]<<
863 "Dist="<<k<<
864 "Dist0="<<distances[0]<<
865 "Dist1="<<distances[1]<<
866 "Dist2="<<distances[2]<<
867 "TDC="<<tdc<<
868 "\n";
869 }
870 return distances[3];
871}
872
873
874
875void AliTOFtrackerMI::GetLikelihood(Float_t dy, Float_t dz, const Double_t *cov, AliTOFtrack * /*track*/, Float_t & py, Float_t &pz)
876{
877 //
878 // get likelihood - track covariance taken
879 // 75 % of gauss with expected sigma
880 // 25 % of gauss with extended sigma
881
882 Double_t kMaxSigmaY = 0.6; // ~ 90% of TRD tracks
883 Double_t kMaxSigmaZ = 1.2; // ~ 90% of TRD tracks
884 Double_t kMeanSigmaY = 0.25; // mean TRD sigma
885 Double_t kMeanSigmaZ = 0.5; // mean TRD sigma
886
887
888 Float_t normwidth, normd, p0,p1;
889 Float_t sigmay = TMath::Max(TMath::Sqrt(cov[0]+kMeanSigmaY*kMeanSigmaY),kMaxSigmaY);
890 Float_t sigmaz = TMath::Max(TMath::Sqrt(cov[2]+kMeanSigmaZ*kMeanSigmaZ),kMaxSigmaZ);
891
892 py=0;
893 pz=0;
894 //
895 // py calculation - 75% admixture of original sigma - 25% tails
896 //
897 normwidth = fDy/sigmay;
898 normd = dy/sigmay;
899 p0 = 0.5*(1+TMath::Erf(normd-normwidth*0.5));
900 p1 = 0.5*(1+TMath::Erf(normd+normwidth*0.5));
901 py+= 0.75*(p1-p0);
902 //
903 normwidth = fDy/(3.*sigmay);
904 normd = dy/(3.*sigmay);
905 p0 = 0.5*(1+TMath::Erf(normd-normwidth*0.5));
906 p1 = 0.5*(1+TMath::Erf(normd+normwidth*0.5));
907 py+= 0.25*(p1-p0);
908 //
909 // pz calculation - 75% admixture of original sigma - 25% tails
910 //
911 normwidth = fDz/sigmaz;
912 normd = dz/sigmaz;
913 p0 = 0.5*(1+TMath::Erf(normd-normwidth*0.5));
914 p1 = 0.5*(1+TMath::Erf(normd+normwidth*0.5));
915 pz+= 0.75*(p1-p0);
916 //
917 normwidth = fDz/(3.*sigmaz);
918 normd = dz/(3.*sigmaz);
919 p0 = 0.5*(1+TMath::Erf(normd-normwidth*0.5));
920 p1 = 0.5*(1+TMath::Erf(normd+normwidth*0.5));
921 pz+= 0.25*(p1-p0);
922}