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 fSeeds(new TObjArray(100)),
67 fClustersESD(new TClonesArray("AliESDTOFCluster")),
68 fHitsESD(new TClonesArray("AliESDTOFHit")),
72 //AliTOFtrackerV2 main Ctor
73 for (Int_t ii=0; ii<kMaxCluster; ii++){
75 fWrittenInPos[ii] = -1;
78 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++)
81 for (Int_t ii=0; ii<4; ii++)
84 //_____________________________________________________________________________
85 AliTOFtrackerV2::~AliTOFtrackerV2() {
90 if(!(AliCDBManager::Instance()->GetCacheFlag())){
100 fClustersESD->Delete();
110 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
111 if(fTimesAr[isp]) delete[] fTimesAr[isp];
112 fTimesAr[isp] = NULL;
116 for (Int_t ii=0; ii<4; ii++)
118 delete [] fTrackPos[ii];
120 //_____________________________________________________________________________
121 void AliTOFtrackerV2::GetPidSettings(AliESDpid *esdPID) {
123 // Sets TOF resolution from RecoParams
126 esdPID->GetTOFResponse().SetTimeResolution(fkRecoParam->GetTimeResolution());
128 AliWarning("fkRecoParam not yet set; cannot set PID settings");
130 //_____________________________________________________________________________
131 Int_t AliTOFtrackerV2::PropagateBack(AliESDEvent * const event) {
133 // Gets seeds from ESD event and Match with TOF Clusters
136 //Update the matched ESD tracks
137 // needed in case of call of TOF info before of the selection of matching and in case of no clusters available at all
142 AliInfo("No TOF recPoints to be matched with reconstructed tracks");
146 // initialize RecoParam for current event
147 AliDebug(1,"Initializing params for TOF");
149 fkRecoParam = AliTOFReconstructor::GetRecoParam(); // instantiate reco param from STEER...
151 if (fkRecoParam == 0x0) {
152 AliFatal("No Reco Param found for TOF!!!");
155 //Initialise some counters
160 Int_t ntrk=event->GetNumberOfTracks();
163 //Load ESD tracks into a local Array of ESD Seeds
164 for (Int_t i=0; i<fNseeds; i++){
165 fSeeds->AddLast(event->GetTrack(i));
166 event->GetTrack(i)->SetESDEvent(event);
168 //Prepare ESD tracks candidates for TOF Matching
171 if (fNseeds==0 || fNseedsTOF==0) {
172 AliInfo("No seeds to try TOF match");
174 fClustersESD->Clear();
179 // clusterize before of matching
180 Clusterize(); // fN might change
182 //Second Step with Looser Matching Criterion
185 // switch array from ALL to filter for ESD (moving from all fClusterESD to filtered AliESDEvent->GetESDTOFClusters())
186 TClonesArray* esdTOFHitArr = event->GetESDTOFHits();
187 TClonesArray* esdTOFClArr = event->GetESDTOFClusters();
189 AliInfo(Form("TOF before the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),fClustersESD->GetEntriesFast()));
191 while(esdTOFHitArr->GetEntriesFast()){ // remove all hits
192 delete esdTOFHitArr->RemoveAt(esdTOFHitArr->GetEntriesFast()-1);
194 for(Int_t i=0;i < fHitsESD->GetEntriesFast();i++){
195 AliESDTOFHit *hitToStored = (AliESDTOFHit *) fHitsESD->At(i);
196 AliESDTOFHit *hitNew = new ( (*esdTOFHitArr)[esdTOFHitArr->GetEntriesFast()] ) AliESDTOFHit(*hitToStored);
197 hitNew->SetESDTOFClusterIndex(hitToStored->GetESDTOFClusterIndex());
200 AliInfo(Form("TOF after the matching: hits = %i and clusters = %i",esdTOFHitArr->GetEntriesFast(),esdTOFClArr->GetEntriesFast()));
203 for (Int_t i=0; i<fNseeds; i++) {
204 AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
205 if ((t->GetStatus()&AliESDtrack::kTOFout)==0)continue;
211 fClustersESD->Clear();
216 //_________________________________________________________________________
217 void AliTOFtrackerV2::CollectESD() {
218 //prepare the set of ESD tracks to be matched to clusters in TOF
226 for (Int_t i=0; i<fNseeds; i++) {
228 AliESDtrack *t =(AliESDtrack*)fSeeds->At(i);
229 if ((t->GetStatus()&AliESDtrack::kTPCout)==0)continue;
232 Float_t x = (Float_t)track.GetX(); //New
235 if ( ( (t->GetStatus()&AliESDtrack::kTRDout)!=0 ) ) {
237 AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
239 // TRD 'good' tracks, already propagated at 371 cm
240 if ( x >= AliTOFGeometry::Rmin() ) {
242 if ( track.PropagateToInnerTOF() ) {
244 AliDebug(1,Form(" TRD propagated track till rho = %fcm."
245 " And then the track has been propagated till rho = %fcm.",
246 x, (Float_t)track.GetX()));
248 track.SetSeedIndex(i);
249 t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
253 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
256 else { // TRD 'good' tracks, propagated rho<371cm
258 if ( track.PropagateToInnerTOF() ) {
260 AliDebug(1,Form(" TRD propagated track till rho = %fcm."
261 " And then the track has been propagated till rho = %fcm.",
262 x, (Float_t)track.GetX()));
264 track.SetSeedIndex(i);
265 t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
269 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
274 else { // Propagate the rest of TPCbp
276 AliDebug(1,Form(" Before propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
278 if ( track.PropagateToInnerTOF() ) {
280 AliDebug(1,Form(" TPC propagated track till rho = %fcm."
281 " And then the track has been propagated till rho = %fcm.",
282 x, (Float_t)track.GetX()));
284 track.SetSeedIndex(i);
285 t->UpdateTrackParams(&track,AliESDtrack::kTOFin);
289 AliDebug(1,Form(" After propagation till inner TOF radius, ESDtrackLength=%f, TOFtrackLength=%f",t->GetIntegratedLength(),track.GetIntegratedLength()));
294 AliInfo(Form("Number of TOF seeds = %d (kTRDout371 = %d, kTRDoutLess371 = %d, !kTRDout = %d)",fNseedsTOF,seedsTOF1,seedsTOF3,seedsTOF2));
298 //_________________________________________________________________________
299 void AliTOFtrackerV2::MatchTracks() {
301 //Match ESD tracks to clusters in TOF
304 // Parameters used/regulating the reconstruction
305 static Float_t detDepth=18.;
306 static Float_t padDepth=0.5;
308 const Float_t kSpeedOfLight= 2.99792458e-2; // speed of light [cm/ps]
310 Float_t dY=AliTOFGeometry::XPad();
311 Float_t dZ=AliTOFGeometry::ZPad();
313 Float_t sensRadius = fkRecoParam->GetSensRadius();
314 Float_t stepSize = fkRecoParam->GetStepSize();
315 Float_t scaleFact = fkRecoParam->GetWindowScaleFact();
316 Float_t dyMax=fkRecoParam->GetWindowSizeMaxY();
317 Float_t dzMax=fkRecoParam->GetWindowSizeMaxZ();
318 Float_t dCut=10.;//fkRecoParam->GetDistanceCut(); // This is to be loaded by OCDB. It should be 10cm always.
319 Double_t maxChi2=fkRecoParam->GetMaxChi2TRD();
320 Bool_t timeWalkCorr = fkRecoParam->GetTimeWalkCorr();
321 AliDebug(1,"++++++++++++++TOF Reconstruction Parameters:++++++++++++++");
322 AliDebug(1,Form("TOF sens radius: %f",sensRadius));
323 AliDebug(1,Form("TOF step size: %f",stepSize));
324 AliDebug(1,Form("TOF Window scale factor: %f",scaleFact));
325 AliDebug(1,Form("TOF Window max dy: %f",dyMax));
326 AliDebug(1,Form("TOF Window max dz: %f",dzMax));
327 AliDebug(1,Form("TOF distance Cut: %f",dCut));
328 AliDebug(1,Form("TOF Max Chi2: %f",maxChi2));
329 AliDebug(1,Form("Time Walk Correction? : %d",timeWalkCorr));
331 //Match ESD tracks to clusters in TOF
332 TClonesArray* TOFClArr = fClustersESD; // use temporary array
333 TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
334 TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
336 if(Int_t(detDepth/stepSize) > fNsteps){ // create array for each step
337 // Get the number of propagation steps
338 fNsteps =(Int_t)(detDepth/stepSize);
340 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
341 if(fTimesAr[isp]) delete[] fTimesAr[isp];
344 for(Int_t isp=0;isp < AliPID::kSPECIESC;isp++){
345 fTimesAr[isp] = new Double_t[fNsteps];
348 for (Int_t ii=0; ii<4; ii++)
350 delete [] fTrackPos[ii];
352 for (Int_t ii=0; ii<4; ii++) fTrackPos[ii] = new Float_t[fNsteps];
357 AliDebug(1,Form(" Number of steps to be done %d",fNsteps));
359 AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
362 const Int_t kNclusterMax = 1000; // related to fN value
363 TGeoHMatrix global[kNclusterMax];
364 Int_t clind[kNclusterMax];
365 Bool_t isClusterMatchable[kNclusterMax]; // true if track and cluster were already matched (set to false below upto nc < kNclusterMax)
367 AliTOFtrack trackTOFin;
370 for (Int_t iseed=0; iseed<fSeeds->GetEntriesFast(); iseed++) {
371 AliESDtrack *t =(AliESDtrack*)fSeeds->At(iseed); // ciao replace with loop on ESD + kTOFin
372 if( (t->GetStatus()&AliESDtrack::kTOFin) == 0 ) continue;
376 for (Int_t ii=0; ii<4; ii++)
377 for (Int_t jj=0; jj<fNsteps; jj++) fTrackPos[ii][jj]=0.;
379 for (Int_t ii=0; ii<kNclusterMax; ii++) clind[ii]=-1;
380 for (Int_t ii=0; ii<kNclusterMax; ii++) global[ii] = 0x0;
381 for (Int_t ii=0; ii<kNclusterMax; ii++) isClusterMatchable[ii] = kFALSE;
383 Double_t timesOr[AliPID::kSPECIESC]; t->GetIntegratedTimes(timesOr,AliPID::kSPECIESC); // in ps
385 // Determine a window around the track
387 trackTOFin.GetExternalParameters(x,par);
389 trackTOFin.GetExternalCovariance(cov);
391 if (cov[0]<0. || cov[2]<0.) {
392 AliWarning(Form("Very strange track (%d)! At least one of its covariance matrix diagonal elements is negative!",iseed));
398 ((5*TMath::Sqrt(TMath::Abs(cov[0])) + 0.5*dY + 2.5*TMath::Abs(par[2]))/sensRadius);
401 (5*TMath::Sqrt(TMath::Abs(cov[2])) + 0.5*dZ + 2.5*TMath::Abs(par[3]));
403 Double_t phi=TMath::ATan2(par[0],x) + trackTOFin.GetAlpha();
404 if (phi<-TMath::Pi())phi+=2*TMath::Pi();
405 if (phi>=TMath::Pi())phi-=2*TMath::Pi();
408 //upper limit on window's size.
409 if (dz> dzMax) dz=dzMax;
410 if (dphi*sensRadius> dyMax) dphi=dyMax/sensRadius;
413 // find the clusters in the window of the track
415 for (Int_t k=FindClusterIndex(z-dz); k<TOFClArr->GetEntriesFast(); k++) {
417 if (nc>=kNclusterMax) {
418 AliWarning("No more matchable clusters can be stored! Please, increase the corresponding vectors size.");
422 AliESDTOFCluster *c=(AliESDTOFCluster *) TOFClArr->At(k);
423 if (c->GetZ() > z+dz) break;
424 if (!c->GetStatus()) {
425 AliDebug(1,"Cluster in channel declared bad!");
426 continue; // skip bad channels as declared in OCDB
430 Double_t dph=TMath::Abs(c->GetPhi()-phi);
431 if (dph>TMath::Pi()) dph-=2.*TMath::Pi();
432 if (TMath::Abs(dph)>dphi) continue;
434 Double_t yc=(c->GetPhi() - trackTOFin.GetAlpha())*c->GetR();
435 Double_t p[2]={yc, c->GetZ()};
436 Double_t cov2[3]= {dY*dY/12., 0., dZ*dZ/12.};
437 if (trackTOFin.AliExternalTrackParam::GetPredictedChi2(p,cov2) > maxChi2)continue;
441 Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(c->GetTOFchannel(),ind);
442 AliTOFGeometry::GetVolumePath(ind,path);
443 gGeoManager->cd(path);
444 global[nc] = *gGeoManager->GetCurrentMatrix();
450 AliDebug(1,Form("No available clusters for the track number %d",iseed));
455 AliDebug(1,Form(" Number of available TOF clusters for the track number %d: %d",iseed,nc));
457 //start fine propagation
459 Int_t nStepsDone = 0;
460 for( Int_t istep=0; istep<fNsteps; istep++){
462 // First of all, propagate the track...
463 Float_t xs = AliTOFGeometry::RinTOF()+istep*stepSize;
464 if (!(trackTOFin.PropagateTo(xs))) break;
466 // ...and then, if necessary, rotate the track
467 Double_t ymax = xs*TMath::Tan(0.5*AliTOFGeometry::GetAlpha());
468 Double_t ysect = trackTOFin.GetY();
470 if (!(trackTOFin.Rotate(AliTOFGeometry::GetAlpha()))) break;
471 } else if (ysect <-ymax) {
472 if (!(trackTOFin.Rotate(-AliTOFGeometry::GetAlpha()))) break;
475 Double_t mom = trackTOFin.P();
478 for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
479 Double_t mass=AliPID::ParticleMass(isp);
480 Double_t momz = mom*AliPID::ParticleCharge(isp);
481 fTimesAr[isp][nStepsDone] = stepSize/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
485 for(Int_t isp=0;isp<AliPID::kSPECIESC;isp++){
486 Double_t mass=AliPID::ParticleMass(isp);
487 Double_t momz = mom*AliPID::ParticleCharge(isp);
488 fTimesAr[isp][nStepsDone] = fTimesAr[isp][nStepsDone-1] + (trackTOFin.GetIntegratedLength()-fTrackPos[3][nStepsDone-1])/kSpeedOfLight*TMath::Sqrt(momz*momz+mass*mass)/momz;
492 // store the running point (Globalrf) - fine propagation
494 Double_t r[3]; trackTOFin.GetXYZ(r);
495 fTrackPos[0][nStepsDone]= (Float_t) r[0];
496 fTrackPos[1][nStepsDone]= (Float_t) r[1];
497 fTrackPos[2][nStepsDone]= (Float_t) r[2];
498 fTrackPos[3][nStepsDone]= trackTOFin.GetIntegratedLength();
501 AliDebug(3,Form(" current step %d (%d) - nStepsDone=%d",istep,fNsteps,nStepsDone));
504 if ( nStepsDone == 0 ) {
505 AliDebug(1,Form(" No track points for track number %d",iseed));
510 AliDebug(3,Form(" Number of steps done for the track number %d: %d",iseed,nStepsDone));
513 for (Int_t i=0; i<nc; i++) isClusterMatchable[i] = kFALSE;
517 Bool_t accept = kFALSE;
518 for (Int_t istep=0; istep<nStepsDone; istep++) {
519 Float_t ctrackPos[3];
520 ctrackPos[0] = fTrackPos[0][istep];
521 ctrackPos[1] = fTrackPos[1][istep];
522 ctrackPos[2] = fTrackPos[2][istep];
524 //now see whether the track matches any of the TOF clusters
526 Float_t dist3d[3]={0.,0.,0.};
529 for (Int_t i=0; i<nc; i++) {
531 AliTOFGeometry::IsInsideThePad((TGeoHMatrix*)(&global[i]),ctrackPos,dist3d);
533 // check multiple hit cases
534 AliESDTOFCluster *cmatched=(AliESDTOFCluster *) TOFClArr->At(clind[i]);
536 if(cmatched->GetNTOFhits() > 1){ // correct residual for mean position of the clusters (w.r.t. the first pad/hit)
537 Float_t zmain = cmatched->GetTOFchannel(0)/48;
538 Float_t xmain = cmatched->GetTOFchannel(0)%48;
539 for(Int_t ihit=1;ihit < cmatched->GetNTOFhits();ihit++){
540 Float_t deltaz = (cmatched->GetTOFchannel(ihit)/48 - zmain) * 3.5;
541 Float_t deltax = (cmatched->GetTOFchannel(ihit)%48 - xmain) * 2.5;
542 dist3d[0] -= deltax / cmatched->GetNTOFhits();
543 dist3d[2] -= deltaz / cmatched->GetNTOFhits();
548 /* if track is inside this cluster set flags which will then
549 * inhibit to add track points for the other clusters */
550 Float_t yLoc = dist3d[1];
551 Float_t rLoc = TMath::Sqrt(dist3d[0]*dist3d[0]+dist3d[2]*dist3d[2]);
552 accept = (TMath::Abs(yLoc)<padDepth*0.5 && rLoc<dCut);
555 /* add point everytime that:
556 * - the tracks is within dCut from the cluster
560 Double_t timesCurrent[AliPID::kSPECIESC];
561 AliDebug(3,Form(" Momentum for track %d -> %f", iseed,t->P()));
562 for (Int_t j=0;j<AliPID::kSPECIESC;j++) {
563 timesCurrent[j] = timesOr[j] + fTimesAr[j][istep];
567 if (TMath::Abs(dist3d[1])<stepSize && !isClusterMatchable[i]) {
568 isClusterMatchable[i] = kTRUE;
570 Int_t currentpos = esdTOFClArr->GetEntriesFast(); // position of cluster in ESD
571 if(fWrittenInPos[clind[i]] != -1){
572 currentpos = fWrittenInPos[clind[i]];
573 cmatched = (AliESDTOFCluster *) esdTOFClArr->At(currentpos); // update the new one in the ESDEvent
575 else{ // add as a new cluster in the ESD TClonesArray
576 AliESDTOFCluster *clnew = new( (*esdTOFClArr)[currentpos] ) AliESDTOFCluster(*cmatched);
577 clnew->SetEvent(fEvent);
578 clnew->SetESDID(currentpos);
580 // remap also TOF hit in the filtered array
581 for(Int_t ii=0;ii < cmatched->GetNTOFhits();ii++){
582 Int_t index = cmatched->GetHitIndex(ii);
583 AliESDTOFHit *hitOld = (AliESDTOFHit *) esdTOFHitArr->At(index);
584 Int_t index_new = fHitsESD->GetEntriesFast();
585 AliESDTOFHit *hitNew = new( (*fHitsESD)[index_new] ) AliESDTOFHit(*hitOld);
586 hitNew->SetESDTOFClusterIndex(currentpos);
587 clnew->SetHitIndex(ii,index_new);
590 fWrittenInPos[clind[i]] = currentpos;
591 cmatched = clnew; // update the new one added to the ESDEvent
594 if(cmatched->GetNMatchableTracks() < AliESDTOFCluster::kMaxMatches){
595 cmatched->Update(t->GetID(),dist3d[0],dist3d[1],dist3d[2],fTrackPos[3][istep],timesCurrent);//x,y,z -> tracking RF
596 t->AddTOFcluster(currentpos);
597 t->SetStatus(AliESDtrack::kTOFout);
600 AliDebug(2,Form(" dist3dLoc[0] = %f, dist3dLoc[1] = %f, dist3dLoc[2] = %f ",dist3d[0],dist3d[1],dist3d[2]));
607 } //end for on the clusters
608 } //end for on the steps
612 AliDebug(1,Form(" No matchable track points for the track number %d",iseed));
617 AliDebug(1,Form(" Number of track points for the track number %d: %d",iseed,nfound));
619 Int_t nMatchedClusters = t->GetNTOFclusters();
621 if (nMatchedClusters==0) {
622 AliDebug(1,Form("Reconstructed track %d doesn't match any TOF cluster", iseed));
627 AliDebug(1,Form(" %d - matched (%d)",iseed,nMatchedClusters));
632 AliTOFcluster cTOF = AliTOFcluster(volIdClus,
633 (Float_t)posClus[0],(Float_t)posClus[1],(Float_t)posClus[2],
634 (Float_t)covClus[0],(Float_t)covClus[1],(Float_t)covClus[2],
635 (Float_t)covClus[3],(Float_t)covClus[4],(Float_t)covClus[5],
636 tofLabels,volIndices,parClu,kTRUE,index[i]);
638 // Fill the track residual histograms.
639 FillResiduals(trackTOFin,c,kFALSE);
644 //_________________________________________________________________________
645 Int_t AliTOFtrackerV2::LoadClusters(TTree *cTree) {
646 //--------------------------------------------------------------------
647 //This function loads the TOF clusters
648 //--------------------------------------------------------------------
650 TBranch *branch=cTree->GetBranch("TOF");
652 AliError("can't get the branch with the TOF clusters !");
656 static TClonesArray dummy("AliTOFcluster",10000);
658 TClonesArray *clusters=&dummy;
659 branch->SetAddress(&clusters);
662 Int_t ncl =clusters->GetEntriesFast();
663 AliInfo(Form("Number of clusters: %d",ncl));
665 fN = ncl; // set cluster counter
667 for(Int_t i=0;i < ncl;i++) // reset position of clusters in ESD
668 fWrittenInPos[i] = -1;
674 for (Int_t i=0; i<ncl; i++) {
675 AliTOFcluster *c=(AliTOFcluster*)clusters->UncheckedAt(i);
679 ind[0]=c->GetDetInd(0);
680 ind[1]=c->GetDetInd(1);
681 ind[2]=c->GetDetInd(2);
682 ind[3]=c->GetDetInd(3);
683 ind[4]=c->GetDetInd(4);
684 Int_t calindex = AliTOFGeometry::GetIndex(ind);
685 Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
693 //_________________________________________________________________________
694 void AliTOFtrackerV2::UnloadClusters() {
695 //--------------------------------------------------------------------
696 //This function unloads TOF clusters
697 //--------------------------------------------------------------------
699 // don't delete TOF clusters here because they should be written
702 //_________________________________________________________________________
703 Int_t AliTOFtrackerV2::FindClusterIndex(Double_t z) const {
704 //--------------------------------------------------------------------
705 // This function returns the index of the nearest cluster
706 //--------------------------------------------------------------------
707 TClonesArray* TOFClArr = fClustersESD;; // use temporary array
708 Int_t n = TOFClArr->GetEntriesFast();
711 if (z <= ((AliESDTOFCluster *) TOFClArr->At(0))->GetZ()) return 0;
712 if (z > ((AliESDTOFCluster *) TOFClArr->At(n-1))->GetZ()) return n;
713 Int_t b=0, e=n-1, m=(b+e)/2;
714 for (; b<e; m=(b+e)/2) {
715 if (z > ((AliESDTOFCluster *) TOFClArr->At(m))->GetZ()) b=m+1;
721 //_________________________________________________________________________
722 Bool_t AliTOFtrackerV2::GetTrackPoint(Int_t index, AliTrackPoint& p) const
724 // Get track space point with index i
725 // Coordinates are in the global system
726 TClonesArray* esdTOFClArr = fEvent->GetESDTOFClusters();
727 AliESDTOFCluster *cl = (AliESDTOFCluster *) esdTOFClArr->At(index);
730 xyz[0] = cl->GetR()*TMath::Cos(cl->GetPhi());
731 xyz[1] = cl->GetR()*TMath::Sin(cl->GetPhi());
733 Float_t phiangle = (Int_t(cl->GetPhi()*TMath::RadToDeg()/20.)+0.5)*20.*TMath::DegToRad();
734 Float_t sinphi = TMath::Sin(phiangle), cosphi = TMath::Cos(phiangle);
735 Int_t tofChannel=cl->GetTOFchannel();
736 Int_t ind[5]; AliTOFGeometry::GetVolumeIndices(tofChannel,ind);
737 Float_t tiltangle = AliTOFGeometry::GetAngles(ind[1],ind[2])*TMath::DegToRad();
738 Float_t sinth = TMath::Sin(tiltangle), costh = TMath::Cos(tiltangle);
739 Float_t sigmay2 = AliTOFGeometry::XPad()*AliTOFGeometry::XPad()/12.;
740 Float_t sigmaz2 = AliTOFGeometry::ZPad()*AliTOFGeometry::ZPad()/12.;
742 cov[0] = sinphi*sinphi*sigmay2 + cosphi*cosphi*sinth*sinth*sigmaz2;
743 cov[1] = -sinphi*cosphi*sigmay2 + sinphi*cosphi*sinth*sinth*sigmaz2;
744 cov[2] = -cosphi*sinth*costh*sigmaz2;
745 cov[3] = cosphi*cosphi*sigmay2 + sinphi*sinphi*sinth*sinth*sigmaz2;
746 cov[4] = -sinphi*sinth*costh*sigmaz2;
747 cov[5] = costh*costh*sigmaz2;
748 p.SetXYZ(xyz[0],xyz[1],xyz[2],cov);
750 // Detector numbering scheme
751 Int_t nSector = AliTOFGeometry::NSectors();
752 Int_t nPlate = AliTOFGeometry::NPlates();
753 Int_t nStripA = AliTOFGeometry::NStripA();
754 Int_t nStripB = AliTOFGeometry::NStripB();
755 Int_t nStripC = AliTOFGeometry::NStripC();
757 Int_t isector = ind[0];//cl->GetDetInd(0);
758 if (isector >= nSector)
759 AliError(Form("Wrong sector number in TOF (%d) !",isector));
760 Int_t iplate = ind[1];//cl->GetDetInd(1);
761 if (iplate >= nPlate)
762 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
763 Int_t istrip = ind[2];//cl->GetDetInd(2);
765 Int_t stripOffset = 0;
771 stripOffset = nStripC;
774 stripOffset = nStripC+nStripB;
777 stripOffset = nStripC+nStripB+nStripA;
780 stripOffset = nStripC+nStripB+nStripA+nStripB;
783 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
787 Int_t idet = (2*(nStripC+nStripB)+nStripA)*isector +
790 UShort_t volid = AliGeomManager::LayerToVolUID(AliGeomManager::kTOF,idet);
791 p.SetVolumeID((UShort_t)volid);
794 //_________________________________________________________________________
796 Float_t AliTOFtrackerV2::CorrectTimeWalk( Float_t dist, Float_t tof) const {
798 //dummy, for the moment
800 if(dist<AliTOFGeometry::ZPad()*0.5){
802 //place here the actual correction
808 //_________________________________________________________________________
809 void AliTOFtrackerV2::Clusterize(){
812 // Load 1 hit in 1 cluster (ESD)
813 TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
814 TClonesArray* esdTOFHitArr = fEvent->GetESDTOFHits();
816 if(TOFClArr->GetEntriesFast()) TOFClArr->Clear();
817 if(esdTOFHitArr->GetEntriesFast()) esdTOFHitArr->Clear();
819 AliESDTOFCluster *c1 = NULL;
820 AliESDTOFCluster *c2 = NULL;
822 for(Int_t i=0; i < fN;i++){
823 AliTOFcluster *c = fClusters[i];
825 ind[0]=c->GetDetInd(0);
826 ind[1]=c->GetDetInd(1);
827 ind[2]=c->GetDetInd(2);
828 ind[3]=c->GetDetInd(3);
829 ind[4]=c->GetDetInd(4);
830 Int_t calindex = AliTOFGeometry::GetIndex(ind);
831 Int_t tofLabels[3]={c->GetLabel(0),c->GetLabel(1),c->GetLabel(2)};
833 new ( (*esdTOFHitArr)[i] ) AliESDTOFHit( AliTOFGeometry::TdcBinWidth()*c->GetTDC(),
834 AliTOFGeometry::TdcBinWidth()*c->GetTDCRAW(),
835 AliTOFGeometry::ToTBinWidth()*c->GetToT()*1E-3,
836 calindex,tofLabels,c->GetL0L1Latency(),
837 c->GetDeltaBC(),i,c->GetZ(),c->GetR(),c->GetPhi() );
839 c1 = new( (*TOFClArr)[i] ) AliESDTOFCluster(i);
840 c1->SetEvent(fEvent);
841 c1->SetStatus( c->GetStatus() );
844 // register hits in the cluster
845 c1->AddESDTOFHitIndex(i);
848 // start to merge clusters
849 Int_t chan1,chan2,chan3;
851 Int_t iphi,iphi2,iphi3;
852 Int_t ieta,ieta2,ieta3;
854 for(Int_t i=0; i < TOFClArr->GetEntriesFast()-1;i++){
855 c1=(AliESDTOFCluster *) TOFClArr->At(i);
856 if(!c1->GetStatus()) continue;
858 chan1 = c1->GetTOFchannel(0);
859 AliTOFGeometry::GetVolumeIndices(chan1, detId); // Get volume index from channel index
861 ieta = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
862 if(detId[1]/*module*/ == 0) ieta += 0;
863 else if(detId[1] == 1) ieta += 38;
864 else if(detId[1] == 2) ieta += 76;
865 else if(detId[1] == 3) ieta += 106;
866 else if(detId[1] == 4) ieta += 144;
867 iphi = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
870 if(c1->GetNTOFhits() > 1){
871 chan2 = c1->GetTOFchannel(1);
872 AliTOFGeometry::GetVolumeIndices(chan2, detId); // Get volume index from channel index
874 ieta2 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
875 if(detId[1]/*module*/ == 0) ieta2 += 0;
876 else if(detId[1] == 1) ieta2 += 38;
877 else if(detId[1] == 2) ieta2 += 76;
878 else if(detId[1] == 3) ieta2 += 106;
879 else if(detId[1] == 4) ieta2 += 144;
880 iphi2 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
887 // 1 and 2 belong now to the first cluster, 3 to the second one
890 for(Int_t j=i+1; j < TOFClArr->GetEntriesFast();j++){
891 c2=(AliESDTOFCluster *) TOFClArr->At(j);
892 if(!c2->GetStatus()) continue;
894 chan3 = c2->GetTOFchannel();
896 // check if the two TOF hits are in the same strip
898 if(strip1 != strip2) continue;
900 AliTOFGeometry::GetVolumeIndices(chan3, detId); // Get volume index from channel index
901 ieta3 = detId[2]/*strip*/*2 + detId[3]/*pad Z*/;
902 if(detId[1]/*module*/ == 0) ieta3 += 0;
903 else if(detId[1] == 1) ieta3 += 38;
904 else if(detId[1] == 2) ieta3 += 76;
905 else if(detId[1] == 3) ieta3 += 106;
906 else if(detId[1] == 4) ieta3 += 144;
907 iphi3 = detId[0]/*phi sector*/*48 + detId[4]/*pad x*/;
910 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")
912 // check if the fired pad are close in space
913 if((TMath::Abs(iphi-iphi3)>1 && TMath::Abs(iphi2-iphi3)>1) || (TMath::Abs(ieta-ieta3)>1 && TMath::Abs(ieta2-ieta3)>1))
914 continue; // double checks
916 // check if the TOF time are close enough to be merged
917 if(TMath::Abs(c1->GetTime() - c2->GetTime()) > 500/*in ps*/) continue;
922 // new hit is added as a second hit for the first cluster
929 void AliTOFtrackerV2::MergeClusters(Int_t i,Int_t j){
930 TClonesArray* TOFClArr = fClustersESD;// use a temporary copy //fEvent->GetESDTOFClusters();
933 AliInfo("No TOF cluster mergine possible (cannot merge a cluster with itself)");
937 if(i > j){ // check right order
943 Int_t last = TOFClArr->GetEntriesFast()-1;
946 AliInfo("No TOF cluster mergine possible (cluster not available)");
950 AliESDTOFCluster *c1 = (AliESDTOFCluster *) TOFClArr->At(i);
951 AliESDTOFCluster *c2 = (AliESDTOFCluster *) TOFClArr->At(j);
953 if(c2->GetNMatchableTracks()){
954 AliInfo("No TOF cluster mergine possible (cluster already matched)");
955 return; // cannot merge a cluster already matched
958 Int_t nhit1 = c1->GetNTOFhits();
959 Int_t nhit2 = c2->GetNTOFhits();
961 if(nhit1+nhit2 >= AliESDTOFCluster::kMaxHits)
963 AliInfo("No TOF cluster mergine possible (too many hits)");
967 for(Int_t k=0;k < nhit2;k++){// add hits in c2 to c1
968 c1->AddESDTOFHitIndex(c2->GetHitIndex(k));
970 // ID re-setting for hits not needed (done when matching is found)
973 // remove c2 from array
974 if(j == last) delete TOFClArr->RemoveAt(j);
976 for(Int_t ii=j;ii < last;ii++){
977 AliESDTOFCluster *old= (AliESDTOFCluster *) TOFClArr->At(ii);
978 if (!old) {AliFatal(Form("NULL pointer for TOF cluster %d",ii));}
979 AliESDTOFCluster *replace= (AliESDTOFCluster *) TOFClArr->At(ii+1);
980 if (!replace) {AliFatal(Form("NULL pointer for TOF cluster %d",ii+1));}
984 delete TOFClArr->RemoveAt(last);