1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //--------------------------------------------------------------------//
18 // AliTOFtrackerV2 Class //
19 // Task: Perform association of the ESD tracks to TOF Clusters //
20 // and Update ESD track with associated TOF Cluster parameters //
22 // -- Authors : S. Arcelli, C. Zampolli (Bologna University and INFN) //
23 // -- Contacts: Annalisa.De.Caro@cern.ch //
24 // -- : Chiara.Zampolli@bo.infn.it //
25 // -- : Silvia.Arcelli@bo.infn.it //
27 //--------------------------------------------------------------------//
32 #include <TClonesArray.h>
33 #include <TObjArray.h>
34 #include <TGeoManager.h>
37 #include "AliGeomManager.h"
38 #include "AliESDtrack.h"
39 #include "AliESDEvent.h"
40 #include "AliESDpid.h"
41 #include "AliESDTOFCluster.h"
43 #include "AliTrackPointArray.h"
44 #include "AliCDBManager.h"
46 #include "AliTOFRecoParam.h"
47 #include "AliTOFReconstructor.h"
48 #include "AliTOFGeometry.h"
49 #include "AliTOFtrackerV2.h"
50 #include "AliTOFtrack.h"
52 #include "AliESDTOFHit.h"
54 extern TGeoManager *gGeoManager;
56 ClassImp(AliTOFtrackerV2)
58 //_____________________________________________________________________________
59 AliTOFtrackerV2::AliTOFtrackerV2():
66 fTracks(new TClonesArray("AliTOFtrack")),
67 fSeeds(new TObjArray(100)),
68 fClustersESD(new TClonesArray("AliESDTOFCluster")),
69 fHitsESD(new TClonesArray("AliESDTOFHit")),
72 //AliTOFtrackerV2 main Ctor
73 for (Int_t ii=0; ii<kMaxCluster; ii++){
75 fWrittenInPos[ii] = -1;
79 //_____________________________________________________________________________
80 AliTOFtrackerV2::~AliTOFtrackerV2() {
85 if(!(AliCDBManager::Instance()->GetCacheFlag())){
100 fClustersESD->Delete();
111 //_____________________________________________________________________________
112 void AliTOFtrackerV2::GetPidSettings(AliESDpid *esdPID) {
114 // Sets TOF resolution from RecoParams
117 esdPID->GetTOFResponse().SetTimeResolution(fkRecoParam->GetTimeResolution());
119 AliWarning("fkRecoParam not yet set; cannot set PID settings");
121 //_____________________________________________________________________________
122 Int_t AliTOFtrackerV2::PropagateBack(AliESDEvent * const event) {
124 // Gets seeds from ESD event and Match with TOF Clusters
127 //Update the matched ESD tracks
128 // needed in case of call of TOF info before of the selection of matching and in case of no clusters available at all
133 AliInfo("No TOF recPoints to be matched with reconstructed tracks");
137 // initialize RecoParam for current event
138 AliDebug(1,"Initializing params for TOF");
140 fkRecoParam = AliTOFReconstructor::GetRecoParam(); // instantiate reco param from STEER...
142 if (fkRecoParam == 0x0) {
143 AliFatal("No Reco Param found for TOF!!!");
146 //Initialise some counters
151 Int_t ntrk=event->GetNumberOfTracks();
154 //Load ESD tracks into a local Array of ESD Seeds
155 for (Int_t i=0; i<fNseeds; i++){
156 fSeeds->AddLast(event->GetTrack(i));
157 event->GetTrack(i)->SetESDEvent(event);
159 //Prepare ESD tracks candidates for TOF Matching
162 if (fNseeds==0 || fNseedsTOF==0) {
163 AliInfo("No seeds to try TOF match");
166 fClustersESD->Clear();
171 // clusterize before of matching
172 Clusterize(); // fN might change
174 //Second Step with Looser Matching Criterion
177 // switch array from ALL to filter for ESD (moving from all fClusterESD to filtered AliESDEvent->GetESDTOFClusters())
178 TClonesArray* esdTOFHitArr = event->GetESDTOFHits();
179 TClonesArray* esdTOFClArr = event->GetESDTOFClusters();
181 AliInfo(Form("TOF before the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),fClustersESD->GetEntriesFast()));
183 while(esdTOFHitArr->GetEntriesFast()){ // remove all hits
184 delete esdTOFHitArr->RemoveAt(esdTOFHitArr->GetEntriesFast()-1);
186 for(Int_t i=0;i < fHitsESD->GetEntriesFast();i++){
187 AliESDTOFHit *hitToStored = (AliESDTOFHit *) fHitsESD->At(i);
188 AliESDTOFHit *hitNew = new ( (*esdTOFHitArr)[esdTOFHitArr->GetEntriesFast()] ) AliESDTOFHit(*hitToStored);
189 hitNew->SetESDTOFClusterIndex(hitToStored->GetESDTOFClusterIndex());
192 AliInfo(Form("TOF after the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),esdTOFClArr->GetEntriesFast()));
195 for (Int_t i=0; i<fNseeds; i++) {
196 AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
197 if ((t->GetStatus()&AliESDtrack::kTOFout)==0)continue;
204 fClustersESD->Clear();
209 //_________________________________________________________________________
210 void AliTOFtrackerV2::CollectESD() {
211 //prepare the set of ESD tracks to be matched to clusters in TOF
217 TClonesArray &aTOFTrack = *fTracks;
218 for (Int_t i=0; i<fNseeds; i++) {
220 AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
221 if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
223 AliTOFtrack *track = new AliTOFtrack(*t); // New
224 Float_t x = (Float_t)track->GetX(); //New
227 if ( ( (t->GetStatus()&AliESDtrack::kTRDout)!=0 ) ) {
229 AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
231 // TRD 'good' tracks, already propagated at 371 cm
232 if ( x >= AliTOFGeometry::Rmin() ) {
234 if ( track->PropagateToInnerTOF() ) {
236 AliDebug(1,Form(" TRD propagated track till rho = %fcm."
237 " And then the track has been propagated till rho = %fcm.",
238 x, (Float_t)track->GetX()));
240 track->SetSeedIndex(i);
241 t->UpdateTrackParams(track,AliESDtrack::kTOFin);
242 new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
246 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
251 else { // TRD 'good' tracks, propagated rho<371cm
253 if ( track->PropagateToInnerTOF() ) {
255 AliDebug(1,Form(" TRD propagated track till rho = %fcm."
256 " And then the track has been propagated till rho = %fcm.",
257 x, (Float_t)track->GetX()));
259 track->SetSeedIndex(i);
260 t->UpdateTrackParams(track,AliESDtrack::kTOFin);
261 new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
265 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
273 else { // Propagate the rest of TPCbp
275 AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
277 if ( track->PropagateToInnerTOF() ) {
279 AliDebug(1,Form(" TPC propagated track till rho = %fcm."
280 " And then the track has been propagated till rho = %fcm.",
281 x, (Float_t)track->GetX()));
283 track->SetSeedIndex(i);
284 t->UpdateTrackParams(track,AliESDtrack::kTOFin);
285 new(aTOFTrack[fNseedsTOF]) AliTOFtrack(*track);
289 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track->GetIntegratedLength()));
295 AliInfo(Form("Number of TOF seeds = %d (kTRDout371 = %d, kTRDoutLess371 = %d, !kTRDout = %d)",fNseedsTOF,seedsTOF1,seedsTOF3,seedsTOF2));
297 // Sort according uncertainties on track position
302 //_________________________________________________________________________
303 void AliTOFtrackerV2::MatchTracks() {
305 //Match ESD tracks to clusters in TOF
308 // Parameters used/regulating the reconstruction
309 static Float_t detDepth=18.;
310 static Float_t padDepth=0.5;
312 const Float_t kSpeedOfLight= 2.99792458e-2; // speed of light [cm/ps]
314 Float_t dY=AliTOFGeometry::XPad();
315 Float_t dZ=AliTOFGeometry::ZPad();
317 Float_t sensRadius = fkRecoParam->GetSensRadius();
318 Float_t stepSize = fkRecoParam->GetStepSize();
319 Float_t scaleFact = fkRecoParam->GetWindowScaleFact();
320 Float_t dyMax=fkRecoParam->GetWindowSizeMaxY();
321 Float_t dzMax=fkRecoParam->GetWindowSizeMaxZ();
322 Float_t dCut=10.;//fkRecoParam->GetDistanceCut(); // This is to be loaded by OCDB. It should be 10cm always.
323 Double_t maxChi2=fkRecoParam->GetMaxChi2TRD();
324 Bool_t timeWalkCorr = fkRecoParam->GetTimeWalkCorr();
325 AliDebug(1,"++++++++++++++TOF Reconstruction Parameters:++++++++++++++");
326 AliDebug(1,Form("TOF sens radius: %f",sensRadius));
327 AliDebug(1,Form("TOF step size: %f",stepSize));
328 AliDebug(1,Form("TOF Window scale factor: %f",scaleFact));
329 AliDebug(1,Form("TOF Window max dy: %f",dyMax));
330 AliDebug(1,Form("TOF Window max dz: %f",dzMax));
331 AliDebug(1,Form("TOF distance Cut: %f",dCut));
332 AliDebug(1,Form("TOF Max Chi2: %f",maxChi2));
333 AliDebug(1,Form("Time Walk Correction? : %d",timeWalkCorr));
335 //Match ESD tracks to clusters in TOF
336 TClonesArray* TOFClArr = fClustersESD; // use temporary array
337 TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
338 TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
340 // Get the number of propagation steps
341 Int_t nSteps=(Int_t)(detDepth/stepSize);
342 AliDebug(1,Form(" Number of steps to be done %d",nSteps));
344 AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
346 //PH Arrays (moved outside of the loop)
347 Float_t * trackPos[4];
348 for (Int_t ii=0; ii<4; ii++) trackPos[ii] = new Float_t[nSteps];
351 const Int_t kNclusterMax = 1000; // related to fN value
352 TGeoHMatrix global[kNclusterMax];
353 Int_t clind[kNclusterMax];
354 Bool_t isClusterMatchable[kNclusterMax]; // true if track and cluster were already matched (set to false below upto nc < kNclusterMax)
358 for (Int_t iseed=0; iseed<fNseedsTOF; iseed++) {
360 for (Int_t ii=0; ii<4; ii++)
361 for (Int_t jj=0; jj<nSteps; jj++) trackPos[ii][jj]=0.;
363 for (Int_t ii=0; ii<kNclusterMax; ii++) clind[ii]=-1;
364 for (Int_t ii=0; ii<kNclusterMax; ii++) global[ii] = 0x0;
365 for (Int_t ii=0; ii<kNclusterMax; ii++) isClusterMatchable[ii] = kFALSE;
367 AliTOFtrack *track =(AliTOFtrack*)fTracks->UncheckedAt(iseed);
368 AliESDtrack *t =(AliESDtrack*)fSeeds->At(track->GetSeedIndex());
369 AliTOFtrack *trackTOFin = new AliTOFtrack(*track);
371 Double_t timesOr[AliPID::kSPECIESC]; t->GetIntegratedTimes(timesOr,AliPID::kSPECIESC); // in ps
373 // Determine a window around the track
375 trackTOFin->GetExternalParameters(x,par);
377 trackTOFin->GetExternalCovariance(cov);
379 if (cov[0]<0. || cov[2]<0.) {
380 AliWarning(Form("Very strange track (%d)! At least one of its covariance matrix diagonal elements is negative!",iseed));
387 ((5*TMath::Sqrt(TMath::Abs(cov[0])) + 0.5*dY + 2.5*TMath::Abs(par[2]))/sensRadius);
390 (5*TMath::Sqrt(TMath::Abs(cov[2])) + 0.5*dZ + 2.5*TMath::Abs(par[3]));
392 Double_t phi=TMath::ATan2(par[0],x) + trackTOFin->GetAlpha();
393 if (phi<-TMath::Pi())phi+=2*TMath::Pi();
394 if (phi>=TMath::Pi())phi-=2*TMath::Pi();
397 //upper limit on window's size.
398 if (dz> dzMax) dz=dzMax;
399 if (dphi*sensRadius> dyMax) dphi=dyMax/sensRadius;
402 // find the clusters in the window of the track
404 for (Int_t k=FindClusterIndex(z-dz); k<TOFClArr->GetEntriesFast(); k++) {
406 if (nc>=kNclusterMax) {
407 AliWarning("No more matchable clusters can be stored! Please, increase the corresponding vectors size.");
411 AliESDTOFCluster *c=(AliESDTOFCluster *) TOFClArr->At(k);
412 if (c->GetZ() > z+dz) break;
413 if (!c->GetStatus()) {
414 AliDebug(1,"Cluster in channel declared bad!");
415 continue; // skip bad channels as declared in OCDB
419 Double_t dph=TMath::Abs(c->GetPhi()-phi);
420 if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
421 if (TMath::Abs(dph)>dphi) continue;
423 Double_t yc=(c->GetPhi() - trackTOFin->GetAlpha())*c->GetR();
424 Double_t p[2]={yc, c->GetZ()};
425 Double_t cov2[3]= {dY*dY/12., 0., dZ*dZ/12.};
426 if (trackTOFin->AliExternalTrackParam::GetPredictedChi2(p,cov2) > maxChi2)continue;
430 Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(c->GetTOFchannel(),ind);
431 AliTOFGeometry::GetVolumePath(ind,path);
432 gGeoManager->cd(path);
433 global[nc] = *gGeoManager->GetCurrentMatrix();
439 AliDebug(1,Form("No available clusters for the track number %d",iseed));
445 AliDebug(1,Form(" Number of available TOF clusters for the track number %d: %d",iseed,nc));
447 //start fine propagation
449 Double_t *times[AliPID::kSPECIESC];
450 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
451 times[isp] = new Double_t[nSteps];
454 Int_t nStepsDone = 0;
455 for( Int_t istep=0; istep<nSteps; istep++){
457 // First of all, propagate the track...
458 Float_t xs = AliTOFGeometry::RinTOF()+istep*stepSize;
459 if (!(trackTOFin->PropagateTo(xs))) break;
461 // ...and then, if necessary, rotate the track
462 Double_t ymax = xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
463 Double_t ysect = trackTOFin->GetY();
465 if (!(trackTOFin->Rotate(AliTOFGeometry::GetAlpha()))) break;
466 } else if (ysect <-ymax) {
467 if (!(trackTOFin->Rotate(-AliTOFGeometry::GetAlpha()))) break;
470 Double_t mom = trackTOFin->P();
473 for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
474 Double_t mass=AliPID::ParticleMass(isp);
475 Double_t momz = mom*AliPID::ParticleCharge(isp);
476 times[isp][nStepsDone] = stepSize/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
480 for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
481 Double_t mass=AliPID::ParticleMass(isp);
482 Double_t momz = mom*AliPID::ParticleCharge(isp);
483 times[isp][nStepsDone] = times[isp][nStepsDone-1] + (trackTOFin->GetIntegratedLength()-trackPos[3][nStepsDone-1])/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
487 // store the running point (Globalrf) - fine propagation
489 Double_t r[3]; trackTOFin->GetXYZ(r);
490 trackPos[0][nStepsDone]= (Float_t) r[0];
491 trackPos[1][nStepsDone]= (Float_t) r[1];
492 trackPos[2][nStepsDone]= (Float_t) r[2];
493 trackPos[3][nStepsDone]= trackTOFin->GetIntegratedLength();
496 AliDebug(3,Form(" current step %d (%d) - nStepsDone=%d",istep,nSteps,nStepsDone));
499 if ( nStepsDone == 0 ) {
500 AliDebug(1,Form(" No track points for track number %d",iseed));
506 AliDebug(3,Form(" Number of steps done for the track number %d: %d",iseed,nStepsDone));
509 for (Int_t i=0; i<nc; i++) isClusterMatchable[i] = kFALSE;
513 Bool_t accept = kFALSE;
514 for (Int_t istep=0; istep<nStepsDone; istep++) {
515 Float_t ctrackPos[3];
516 ctrackPos[0] = trackPos[0][istep];
517 ctrackPos[1] = trackPos[1][istep];
518 ctrackPos[2] = trackPos[2][istep];
520 //now see whether the track matches any of the TOF clusters
522 Float_t dist3d[3]={0.,0.,0.};
525 for (Int_t i=0; i<nc; i++) {
527 AliTOFGeometry::IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
530 /* if track is inside this cluster set flags which will then
531 * inhibit to add track points for the other clusters */
532 Float_t yLoc = dist3d[1];
533 Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
534 accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
537 /* add point everytime that:
538 * - the tracks is within dCut from the cluster
542 Double_t timesCurrent[AliPID::kSPECIESC];
543 AliDebug(3,Form(" Momentum for track %d -> %f", iseed,t->P()));
544 for (Int_t j=0;j<AliPID::kSPECIESC;j++) {
545 timesCurrent[j] = timesOr[j] + times[j][istep];
549 if (TMath::Abs(dist3d[1])<stepSize && !isClusterMatchable[i]) {
550 isClusterMatchable[i] = kTRUE;
552 AliESDTOFCluster *cmatched=(AliESDTOFCluster *) TOFClArr->At(clind[i]);
554 Int_t currentpos = esdTOFClArr->GetEntriesFast(); // position of cluster in ESD
555 if(fWrittenInPos[clind[i]] != -1){
556 currentpos = fWrittenInPos[clind[i]];
557 cmatched = (AliESDTOFCluster *) esdTOFClArr->At(currentpos); // update the new one in the ESDEvent
559 else{ // add as a new cluster in the ESD TClonesArray
560 AliESDTOFCluster *clnew = new( (*esdTOFClArr)[currentpos] ) AliESDTOFCluster(*cmatched);
561 clnew->SetEvent(fEvent);
562 clnew->SetESDID(currentpos);
564 // remap also TOF hit in the filtered array
565 for(Int_t ii=0;ii < cmatched->GetNTOFhits();ii++){
566 Int_t index = cmatched->GetHitIndex(ii);
567 AliESDTOFHit *hitOld = (AliESDTOFHit *) esdTOFHitArr->At(index);
568 Int_t index_new = fHitsESD->GetEntriesFast();
569 AliESDTOFHit *hitNew = new( (*fHitsESD)[index_new] ) AliESDTOFHit(*hitOld);
570 hitNew->SetESDTOFClusterIndex(currentpos);
571 clnew->SetHitIndex(ii,index_new);
574 fWrittenInPos[clind[i]] = currentpos;
575 cmatched = clnew; // update the new one added to the ESDEvent
578 if(cmatched->GetNMatchableTracks() < AliESDTOFCluster::kMaxMatches){
579 printf("dist3d = %f,%f,%f\n",dist3d[0],dist3d[1],dist3d[2]);
580 cmatched->Update(t->GetID(),dist3d[0],dist3d[1],dist3d[2],trackPos[3][istep],timesCurrent);//x,y,z -> tracking RF
581 t->AddTOFcluster(currentpos);
582 t->SetStatus(AliESDtrack::kTOFout);
585 AliDebug(2,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
592 } //end for on the clusters
593 } //end for on the steps
595 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
601 AliDebug(1,Form(" No matchable track points for the track number %d",iseed));
607 AliDebug(1,Form(" Number of track points for the track number %d: %d",iseed,nfound));
609 Int_t nMatchedClusters = t->GetNTOFclusters();
611 if (nMatchedClusters==0) {
612 AliDebug(1,Form("Reconstructed track %d doesn't match any TOF cluster", iseed));
618 AliDebug(1,Form(" %d - matched (%d)",track->GetSeedIndex()/*iseed*/,nMatchedClusters));
623 AliTOFcluster cTOF = AliTOFcluster(volIdClus,
624 (Float_t)posClus[0],(Float_t)posClus[1],(Float_t)posClus[2],
625 (Float_t)covClus[0],(Float_t)covClus[1],(Float_t)covClus[2],
626 (Float_t)covClus[3],(Float_t)covClus[4],(Float_t)covClus[5],
627 tofLabels,volIndices,parClu,kTRUE,index[i]);
629 // Fill the track residual histograms.
630 FillResiduals(trackTOFin,c,kFALSE);
637 for (Int_t ii=0; ii<4; ii++) delete [] trackPos[ii];
640 //_________________________________________________________________________
641 Int_t AliTOFtrackerV2::LoadClusters(TTree *cTree) {
642 //--------------------------------------------------------------------
643 //This function loads the TOF clusters
644 //--------------------------------------------------------------------
646 TBranch *branch=cTree->GetBranch("TOF");
648 AliError("can't get the branch with the TOF clusters !");
652 static TClonesArray dummy("AliTOFcluster",10000);
654 TClonesArray *clusters=&dummy;
655 branch->SetAddress(&clusters);
658 Int_t ncl =clusters->GetEntriesFast();
659 AliInfo(Form("Number of clusters: %d",ncl));
661 fN = ncl; // set cluster counter
663 for(Int_t i=0;i < ncl;i++) // reset position of clusters in ESD
664 fWrittenInPos[i] = -1;
670 for (Int_t i=0; i<ncl; i++) {
671 AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
675 ind[0]=c->GetDetInd(0);
676 ind[1]=c->GetDetInd(1);
677 ind[2]=c->GetDetInd(2);
678 ind[3]=c->GetDetInd(3);
679 ind[4]=c->GetDetInd(4);
680 Int_t calindex = AliTOFGeometry::GetIndex(ind);
681 Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
689 //_________________________________________________________________________
690 void AliTOFtrackerV2::UnloadClusters() {
691 //--------------------------------------------------------------------
692 //This function unloads TOF clusters
693 //--------------------------------------------------------------------
695 // don't delete TOF clusters here because they should be written
698 //_________________________________________________________________________
699 Int_t AliTOFtrackerV2::FindClusterIndex(Double_t z) const {
700 //--------------------------------------------------------------------
701 // This function returns the index of the nearest cluster
702 //--------------------------------------------------------------------
703 TClonesArray* TOFClArr = fClustersESD;; // use temporary array
704 Int_t n = TOFClArr->GetEntriesFast();
707 if (z <= ((AliESDTOFCluster *) TOFClArr->At(0))->GetZ()) return 0;
708 if (z > ((AliESDTOFCluster *) TOFClArr->At(n-1))->GetZ()) return n;
709 Int_t b=0, e=n-1, m=(b+e)/2;
710 for (; b<e; m=(b+e)/2) {
711 if (z > ((AliESDTOFCluster *) TOFClArr->At(m))->GetZ()) b=m+1;
717 //_________________________________________________________________________
718 Bool_t AliTOFtrackerV2::GetTrackPoint(Int_t index, AliTrackPoint& p) const
720 // Get track space point with index i
721 // Coordinates are in the global system
722 TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
723 AliESDTOFCluster *cl = (AliESDTOFCluster *) esdTOFClArr->At(index);
726 xyz[0] = cl->GetR()*TMath::Cos(cl->GetPhi());
727 xyz[1] = cl->GetR()*TMath::Sin(cl->GetPhi());
729 Float_t phiangle = (Int_t(cl->GetPhi()*TMath::RadToDeg()/20.)+0.5)*20.*TMath::DegToRad();
730 Float_t sinphi = TMath::Sin(phiangle), cosphi = TMath::Cos(phiangle);
731 Int_t tofChannel=cl->GetTOFchannel();
732 Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(tofChannel,ind);
733 Float_t tiltangle = AliTOFGeometry::GetAngles(ind[1],ind[2])*TMath::DegToRad();
734 Float_t sinth = TMath::Sin(tiltangle), costh = TMath::Cos(tiltangle);
735 Float_t sigmay2 = AliTOFGeometry::XPad()*AliTOFGeometry::XPad()/12.;
736 Float_t sigmaz2 = AliTOFGeometry::ZPad()*AliTOFGeometry::ZPad()/12.;
738 cov[0] = sinphi*sinphi*sigmay2 + cosphi*cosphi*sinth*sinth*sigmaz2;
739 cov[1] = -sinphi*cosphi*sigmay2 + sinphi*cosphi*sinth*sinth*sigmaz2;
740 cov[2] = -cosphi*sinth*costh*sigmaz2;
741 cov[3] = cosphi*cosphi*sigmay2 + sinphi*sinphi*sinth*sinth*sigmaz2;
742 cov[4] = -sinphi*sinth*costh*sigmaz2;
743 cov[5] = costh*costh*sigmaz2;
744 p.SetXYZ(xyz[0],xyz[1],xyz[2],cov);
746 // Detector numbering scheme
747 Int_t nSector = AliTOFGeometry::NSectors();
748 Int_t nPlate = AliTOFGeometry::NPlates();
749 Int_t nStripA = AliTOFGeometry::NStripA();
750 Int_t nStripB = AliTOFGeometry::NStripB();
751 Int_t nStripC = AliTOFGeometry::NStripC();
753 Int_t isector = ind[0];//cl->GetDetInd(0);
754 if (isector >= nSector)
755 AliError(Form("Wrong sector number in TOF (%d) !",isector));
756 Int_t iplate = ind[1];//cl->GetDetInd(1);
757 if (iplate >= nPlate)
758 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
759 Int_t istrip = ind[2];//cl->GetDetInd(2);
761 Int_t stripOffset = 0;
767 stripOffset = nStripC;
770 stripOffset = nStripC+nStripB;
773 stripOffset = nStripC+nStripB+nStripA;
776 stripOffset = nStripC+nStripB+nStripA+nStripB;
779 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
783 Int_t idet = (2*(nStripC+nStripB)+nStripA)*isector +
786 UShort_t volid = AliGeomManager::LayerToVolUID(AliGeomManager::kTOF,idet);
787 p.SetVolumeID((UShort_t)volid);
790 //_________________________________________________________________________
792 Float_t AliTOFtrackerV2::CorrectTimeWalk( Float_t dist, Float_t tof) const {
794 //dummy, for the moment
796 if(dist<AliTOFGeometry::ZPad()*0.5){
798 //place here the actual correction
804 //_________________________________________________________________________
805 void AliTOFtrackerV2::Clusterize(){
808 // Load 1 hit in 1 cluster (ESD)
809 TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
810 TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
812 if(TOFClArr->GetEntriesFast()) TOFClArr->Clear();
813 if(esdTOFHitArr->GetEntriesFast()) esdTOFHitArr->Clear();
815 AliESDTOFCluster *c1 = NULL;
816 AliESDTOFCluster *c2 = NULL;
818 for(Int_t i=0; i < fN;i++){
819 AliTOFcluster *c = fClusters[i];
821 ind[0]=c->GetDetInd(0);
822 ind[1]=c->GetDetInd(1);
823 ind[2]=c->GetDetInd(2);
824 ind[3]=c->GetDetInd(3);
825 ind[4]=c->GetDetInd(4);
826 Int_t calindex = AliTOFGeometry::GetIndex(ind);
827 Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
829 new ( (*esdTOFHitArr)[i] ) AliESDTOFHit( AliTOFGeometry::TdcBinWidth()*c->GetTDC(),
830 AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW(),
831 AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3,
832 calindex,tofLabels,c->GetL0L1Latency(),
833 c->GetDeltaBC(),i,c->GetZ(),c->GetR(),c->GetPhi() );
835 c1 = new( (*TOFClArr)[i] ) AliESDTOFCluster(i);
836 c1->SetEvent(fEvent);
837 c1->SetStatus( c->GetStatus() );
840 // register hits in the cluster
841 c1->AddESDTOFHitIndex(i);
845 // start to merge clusters
846 Int_t chan1,chan2,chan3;
848 Int_t iphi,iphi2,iphi3;
849 Int_t ieta,ieta2,ieta3;
851 for(Int_t i=0; i < TOFClArr->GetEntriesFast()-1;i++){
852 c1=(AliESDTOFCluster *) TOFClArr->At(i);
853 if(!c1->GetStatus()) continue;
855 chan1 = c1->GetTOFchannel();
856 AliTOFGeometry::GetVolumeIndices(chan1, detId); // Get volume index from channel index
858 ieta = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
859 if(detId[1]/*module*/ == 0) ieta += 0;
860 else if(detId[1] == 1) ieta += 38;
861 else if(detId[1] == 2) ieta += 76;
862 else if(detId[1] == 3) ieta += 106;
863 else if(detId[1] == 4) ieta += 144;
864 iphi = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
867 for(Int_t j=i+1; j < TOFClArr->GetEntriesFast();j++){
868 c2=(AliESDTOFCluster *) TOFClArr->At(j);
869 if(!c2->GetStatus()) continue;
871 chan2 = c2->GetTOFchannel();
873 // check if the two TOF hits are in the same strip
875 if(strip1 != strip2) continue;
877 AliTOFGeometry::GetVolumeIndices(chan2, detId); // Get volume index from channel index
878 ieta2 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
879 if(detId[1]/*module*/ == 0) ieta2 += 0;
880 else if(detId[1] == 1) ieta2 += 38;
881 else if(detId[1] == 2) ieta2 += 76;
882 else if(detId[1] == 3) ieta2 += 106;
883 else if(detId[1] == 4) ieta2 += 144;
884 iphi2 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
886 if(ieta2-ieta > 2) j = fN; // because cluster are order along Z so you can skip all the rest, go to the next one ("i+1")
888 // check if the fired pad are close in space
889 if(TMath::Abs(iphi-iphi2)>1 || TMath::Abs(ieta-ieta2)>1) continue;
891 // check if the TOF time are close enough to be merged
892 if(TMath::Abs(c1->GetTime() - c2->GetTime()) > 500/*in ps*/) continue;
897 j = fN; // cluster "i" merged go to the next one ("i+1")
901 // Sort in z after clusterization
902 for(Int_t i=0; i < TOFClArr->GetEntriesFast()-1;i++){
903 c1=(AliESDTOFCluster *) TOFClArr->At(i);
904 for(Int_t j=i+1; j < TOFClArr->GetEntriesFast()-1;j++){
905 c2=(AliESDTOFCluster *) TOFClArr->At(j);
907 if(c2->GetZ() < c1->GetZ()){ // exchange position
908 AliESDTOFCluster ctemp(*c2);
917 // ID re-setting for hits not needed (done when matching is found)
921 // second step of clusterization
922 for(Int_t i=0; i < TOFClArr->GetEntriesFast()-1;i++){
923 c1=(AliESDTOFCluster *) TOFClArr->At(i);
924 if(!c1->GetStatus()) continue;
926 chan1 = c1->GetTOFchannel(0);
927 AliTOFGeometry::GetVolumeIndices(chan1, detId); // Get volume index from channel index
929 ieta = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
930 if(detId[1]/*module*/ == 0) ieta += 0;
931 else if(detId[1] == 1) ieta += 38;
932 else if(detId[1] == 2) ieta += 76;
933 else if(detId[1] == 3) ieta += 106;
934 else if(detId[1] == 4) ieta += 144;
935 iphi = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
938 if(c1->GetNTOFhits() > 1){
939 chan2 = c1->GetTOFchannel(1);
940 AliTOFGeometry::GetVolumeIndices(chan2, detId); // Get volume index from channel index
942 ieta2 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
943 if(detId[1]/*module*/ == 0) ieta2 += 0;
944 else if(detId[1] == 1) ieta2 += 38;
945 else if(detId[1] == 2) ieta2 += 76;
946 else if(detId[1] == 3) ieta2 += 106;
947 else if(detId[1] == 4) ieta2 += 144;
948 iphi2 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
955 // 1 and 2 belong now to the first cluster, 3 to the second one
958 for(Int_t j=i+1; j < TOFClArr->GetEntriesFast();j++){
959 c2=(AliESDTOFCluster *) TOFClArr->At(j);
960 if(!c2->GetStatus()) continue;
962 chan3 = c2->GetTOFchannel();
964 // check if the two TOF hits are in the same strip
966 if(strip1 != strip2) continue;
968 AliTOFGeometry::GetVolumeIndices(chan3, detId); // Get volume index from channel index
969 ieta3 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
970 if(detId[1]/*module*/ == 0) ieta3 += 0;
971 else if(detId[1] == 1) ieta3 += 38;
972 else if(detId[1] == 2) ieta3 += 76;
973 else if(detId[1] == 3) ieta3 += 106;
974 else if(detId[1] == 4) ieta3 += 144;
975 iphi3 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
978 if(ieta3-ieta > 2) j = fN; // because cluster are order along Z so you can skip all the rest, go to the next one ("i+1")
980 // check if the fired pad are close in space
981 if((TMath::Abs(iphi-iphi3)>1 && TMath::Abs(iphi2-iphi3)>1) || (TMath::Abs(ieta-ieta3)>1 && TMath::Abs(ieta2-ieta3)>1))
982 continue; // double checks
984 // check if the TOF time are close enough to be merged
985 if(TMath::Abs(c1->GetTime() - c2->GetTime()) > 500/*in ps*/) continue;
990 j = fN; // cluster "i" merged go to the next one ("i+1")
994 // Sort in z after clusterization
995 for(Int_t i=0; i < TOFClArr->GetEntriesFast()-1;i++){
996 c1=(AliESDTOFCluster *) TOFClArr->At(i);
997 for(Int_t j=i+1; j < TOFClArr->GetEntriesFast()-1;j++){
998 c2=(AliESDTOFCluster *) TOFClArr->At(j);
1000 if(c2->GetZ() < c1->GetZ()){ // exchange position
1001 AliESDTOFCluster ctemp(*c2);
1010 // ID re-setting for hits not needed (done when matching is found)
1015 void AliTOFtrackerV2::MergeClusters(Int_t i,Int_t j){
1016 TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
1017 // TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
1020 AliInfo("No TOF cluster mergine possible (cannot merge a cluster with itself)");
1024 if(i > j){ // check right order
1030 Int_t last = TOFClArr->GetEntriesFast()-1;
1033 AliInfo("No TOF cluster mergine possible (cluster not available)");
1037 AliESDTOFCluster *c1 = (AliESDTOFCluster *) TOFClArr->At(i);
1038 AliESDTOFCluster *c2 = (AliESDTOFCluster *) TOFClArr->At(j);
1040 if(c2->GetNMatchableTracks()){
1041 AliInfo("No TOF cluster mergine possible (cluster already matched)");
1042 return; // cannot merge a cluster already matched
1045 Int_t nhit1 = c1->GetNTOFhits();
1046 Int_t nhit2 = c2->GetNTOFhits();
1048 if(nhit1+nhit2 >= AliESDTOFCluster::kMaxHits)
1050 AliInfo("No TOF cluster mergine possible (too many hits)");
1054 for(Int_t k=0;k < nhit2;k++){// add hits in c2 to c1
1055 c1->AddESDTOFHitIndex(c2->GetHitIndex(k));
1057 // ID re-setting for hits not needed (done when matching is found)
1060 // remove c2 from array
1061 delete TOFClArr->RemoveAt(j);
1063 AliESDTOFCluster *replace= (AliESDTOFCluster *) TOFClArr->At(last);
1064 if (!replace) {AliFatal(Form("NULL pointer for TOF cluster %d",last));}
1065 replace->SetESDID(j);
1066 new ( (*TOFClArr)[j] ) AliESDTOFCluster(*replace);
1067 delete TOFClArr->RemoveAt(last);