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