db8db279b10338117941513627ef417d14195f10
[u/mrichter/AliRoot.git] / TOF / AliTOFClusterFinder.cxx
1 /***************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* 
17
18 Revision 0.03  2005/07/28 A. De Caro
19          Implement public method
20          Raw2Digits(Int_t, AliRawReader *)
21          to convert digits from raw data in MC digits
22          (temporary solution)
23
24 Revision 0.02  2005/07/27 A. De Caro
25          Implement public method
26          Digits2RecPoint(Int_t)
27          to convert digits in clusters
28
29 Revision 0.02  2005/07/26 A. De Caro
30          Implement private methods
31          InsertCluster(AliTOFcluster *)
32          FindClusterIndex(Double_t)
33          originally implemented in AliTOFtracker
34          by S. Arcelli and C. Zampolli
35
36 Revision 0.01  2005/07/25 A. De Caro
37          Implement public methods
38          Digits2RecPoint(AliRawReader *, TTree *)
39          Digits2RecPoint(Int_t, AliRawReader *)
40          to convert raw data in clusters
41  */
42
43 ////////////////////////////////////////////////////////////////
44 //                                                            //
45 //         Class for TOF cluster finder                       //
46 //                                                            //
47 // Starting from Raw Data, create rec points,                 //
48 //                         fill TreeR for TOF,                //
49 //                         write TOF.RecPoints.root file      //
50 //                                                            //
51 ////////////////////////////////////////////////////////////////
52
53 #include "TClonesArray.h"
54 #include "TFile.h"
55 #include "TTree.h"
56
57 #include "AliDAQ.h"
58 #include "AliLoader.h"
59 #include "AliLog.h"
60 #include "AliRawReader.h"
61 #include "AliRunLoader.h"
62
63 #include "AliTOFcalib.h"
64 #include "AliTOFCal.h"
65 #include "AliTOFChannel.h"
66 #include "AliTOFClusterFinder.h"
67 #include "AliTOFcluster.h"
68 #include "AliTOFdigit.h"
69 #include "AliTOFGeometryV5.h"
70 #include "AliTOFGeometry.h"
71 #include "AliTOFRawStream.h"
72
73 extern TFile *gFile;
74
75 ClassImp(AliTOFClusterFinder)
76
77 AliTOFClusterFinder::AliTOFClusterFinder():
78   fRunLoader(0),
79   fTOFLoader(0),
80   fTreeD(0),
81   fTreeR(0),
82   fDigits(new TClonesArray("AliTOFdigit", 4000)),
83   fRecPoints(new TClonesArray("AliTOFcluster", 4000)),
84   fNumberOfTofClusters(0)
85 {
86 //
87 // Constructor
88 //
89
90   fTOFGeometry = new AliTOFGeometryV5();
91   AliInfo("V5 TOF Geometry is taken as the default");
92
93 }
94 //______________________________________________________________________________
95
96 AliTOFClusterFinder::AliTOFClusterFinder(AliRunLoader* runLoader):
97   fRunLoader(runLoader),
98   fTOFLoader(runLoader->GetLoader("TOFLoader")),
99   fTreeD(0),
100   fTreeR(0),
101   fDigits(new TClonesArray("AliTOFdigit", 4000)),
102   fRecPoints(new TClonesArray("AliTOFcluster", 4000)),
103   fNumberOfTofClusters(0)
104 {
105 //
106 // Constructor
107 //
108
109   runLoader->CdGAFile();
110   TFile *in=(TFile*)gFile;
111   in->cd();
112   fTOFGeometry = (AliTOFGeometry*)in->Get("TOFgeometry");
113
114 }
115
116 //------------------------------------------------------------------------
117 AliTOFClusterFinder::AliTOFClusterFinder(const AliTOFClusterFinder &source)
118   :TObject()
119 {
120   // copy constructor
121   this->fDigits=source.fDigits;
122   this->fRecPoints=source.fRecPoints;
123   this->fTOFGeometry=source.fTOFGeometry;
124
125 }
126
127 //------------------------------------------------------------------------
128 AliTOFClusterFinder& AliTOFClusterFinder::operator=(const AliTOFClusterFinder &source)
129 {
130   // ass. op.
131   this->fDigits=source.fDigits;
132   this->fRecPoints=source.fRecPoints;
133   this->fTOFGeometry=source.fTOFGeometry;
134   return *this;
135
136 }
137 //______________________________________________________________________________
138
139 AliTOFClusterFinder::~AliTOFClusterFinder()
140 {
141
142   //
143   // Destructor
144   //
145
146   if (fDigits)
147     {
148       fDigits->Delete();
149       delete fDigits;
150       fDigits=0;
151     }
152   if (fRecPoints)
153     {
154       fRecPoints->Delete();
155       delete fRecPoints;
156       fRecPoints=0;
157     }
158
159   delete fTOFGeometry;
160
161 }
162 //______________________________________________________________________________
163
164 void AliTOFClusterFinder::Digits2RecPoints(Int_t iEvent)
165 {
166   //
167   // Converts digits to recpoints for TOF
168   //
169
170   fRunLoader->GetEvent(iEvent);
171
172   fTreeD = fTOFLoader->TreeD();
173   if (fTreeD == 0x0)
174     {
175       AliFatal("AliTOFClusterFinder: Can not get TreeD");
176     }
177
178   TBranch *branch = fTreeD->GetBranch("TOF");
179   if (!branch) { 
180     AliError("can't get the branch with the TOF digits !");
181     return;
182   }
183
184   TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy;
185   branch->SetAddress(&digits);
186
187   ResetRecpoint();
188
189   fTreeR = fTOFLoader->TreeR();
190   if (fTreeR == 0x0)
191     {
192       fTOFLoader->MakeTree("R");
193       fTreeR = fTOFLoader->TreeR();
194     }
195
196   Int_t bufsize = 32000;
197   fTreeR->Branch("TOF", &fRecPoints, bufsize);
198
199   fTreeD->GetEvent(0);
200   Int_t nDigits = digits->GetEntriesFast();
201   AliDebug(2,Form("Number of TOF digits: %d",nDigits));
202
203   Int_t ii, jj;
204   Int_t dig[5];
205   Float_t g[3];
206   Double_t h[5];
207   Float_t tToT;
208   Double_t tTdcND;
209   for (ii=0; ii<nDigits; ii++) {
210     AliTOFdigit *d = (AliTOFdigit*)digits->UncheckedAt(ii);
211     dig[0]=d->GetSector();
212     dig[1]=d->GetPlate();
213     dig[2]=d->GetStrip();
214     dig[3]=d->GetPadz();
215     dig[4]=d->GetPadx();
216
217     //AliInfo(Form(" %2i  %1i  %2i  %1i  %2i ",dig[0],dig[1],dig[2],dig[3],dig[4]));
218
219     for (jj=0; jj<3; jj++) g[jj] = 0.;
220     fTOFGeometry->GetPos(dig,g);
221
222     h[0] = TMath::Sqrt(g[0]*g[0]+g[1]*g[1]);
223     h[1] = TMath::ATan2(g[1],g[0]);
224     h[2] = g[2];
225     h[3] = d->GetTdc();
226     h[4] = d->GetAdc();
227     tToT = d->GetToT();
228     tTdcND = d->GetTdcND();
229
230     AliTOFcluster *tofCluster = new AliTOFcluster(h,d->GetTracks(),dig,ii,tToT, tTdcND);
231     InsertCluster(tofCluster);
232
233   }
234
235   AliInfo(Form("Number of found clusters: %i", fNumberOfTofClusters));
236
237   CalibrateRecPoint();
238   FillRecPoint();
239
240   fTreeR->Fill();
241   ResetRecpoint();
242
243   fTOFLoader = fRunLoader->GetLoader("TOFLoader");  
244   fTOFLoader->WriteRecPoints("OVERWRITE");
245
246 }
247 //______________________________________________________________________________
248
249 void AliTOFClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
250                                            TTree *clustersTree)
251 {
252   //
253   // Converts RAW data to recpoints for TOF
254   //
255
256   //const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
257   const Int_t kDDL = AliDAQ::NumberOfDdls("TOF");
258
259   ResetRecpoint();
260
261   Int_t bufsize = 32000;
262   clustersTree->Branch("TOF", &fRecPoints, bufsize);
263
264   Int_t ii = 0;
265   Int_t indexDDL = 0;
266
267   Int_t detectorIndex[5];
268   Float_t position[3];
269   Double_t cylindricalPosition[5];
270   Float_t tToT;
271   Double_t tTdcND;
272
273   for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
274
275     rawReader->Reset();
276     AliTOFRawStream tofInput(rawReader);
277     rawReader->Select("TOF", indexDDL, indexDDL);
278
279     while(tofInput.Next()) {
280
281       for (ii=0; ii<5; ii++) detectorIndex[ii] = -1;
282
283       detectorIndex[0] = tofInput.GetSector();
284       detectorIndex[1] = tofInput.GetPlate();
285       detectorIndex[2] = tofInput.GetStrip();
286       detectorIndex[3] = tofInput.GetPadZ();
287       detectorIndex[4] = tofInput.GetPadX();
288       
289       //AliInfo(Form("  %2i  %1i  %2i  %1i  %2i ",detectorIndex[0],detectorIndex[1],detectorIndex[2],detectorIndex[3],detectorIndex[4]));
290
291       if (detectorIndex[0]==-1 ||
292           detectorIndex[1]==-1 ||
293           detectorIndex[2]==-1 ||
294           detectorIndex[3]==-1 ||
295           detectorIndex[4]==-1) continue;
296
297       for (ii=0; ii<3; ii++) position[ii] =  0.;
298
299       fTOFGeometry->GetPos(detectorIndex, position);
300
301       cylindricalPosition[0] = TMath::Sqrt(position[0]*position[0] + position[1]*position[1]);
302       cylindricalPosition[1] = TMath::ATan2(position[1], position[0]);
303       cylindricalPosition[2] = position[2];
304       cylindricalPosition[3] = tofInput.GetTofBin();
305       cylindricalPosition[4] = tofInput.GetToTbin();
306       tToT = tofInput.GetToTbin();
307       tTdcND = -1.;
308       AliTOFcluster *tofCluster = new AliTOFcluster(cylindricalPosition, detectorIndex);
309       tofCluster->SetToT(tToT);
310       tofCluster->SetTDCND(tTdcND);
311       InsertCluster(tofCluster);
312
313     } // while loop
314
315   } // loop on DDL files
316
317   AliInfo(Form("Number of found clusters: %i", fNumberOfTofClusters));
318
319   CalibrateRecPoint();
320   FillRecPoint();
321
322   clustersTree->Fill();
323
324   ResetRecpoint();
325
326 }
327 //______________________________________________________________________________
328
329 void AliTOFClusterFinder::Digits2RecPoints(Int_t iEvent, AliRawReader *rawReader)
330 {
331   //
332   // Converts RAW data to recpoints for TOF
333   //
334
335   //const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
336   const Int_t kDDL = AliDAQ::NumberOfDdls("TOF");
337
338   fRunLoader->GetEvent(iEvent);
339
340   AliDebug(2,Form(" Event number %2i ", iEvent));
341
342   fTreeR = fTOFLoader->TreeR();
343
344   if (fTreeR == 0x0){
345     fTOFLoader->MakeTree("R");
346     fTreeR = fTOFLoader->TreeR();
347   }
348
349   Int_t bufsize = 32000;
350   fTreeR->Branch("TOF", &fRecPoints, bufsize);
351
352   Int_t ii = 0;
353   Int_t indexDDL = 0;
354
355   Int_t detectorIndex[5] = {-1, -1, -1, -1, -1};
356   Float_t position[3];
357   Double_t cylindricalPosition[5];
358   Float_t tToT;
359   Double_t tTdcND;
360
361   for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
362
363     rawReader->Reset();
364     AliTOFRawStream tofInput(rawReader);
365     rawReader->Select("TOF", indexDDL, indexDDL);
366
367     while(tofInput.Next()) {
368
369       for (ii=0; ii<5; ii++) detectorIndex[ii] = -1;
370
371       detectorIndex[0] = (Int_t)tofInput.GetSector();
372       detectorIndex[1] = (Int_t)tofInput.GetPlate();
373       detectorIndex[2] = (Int_t)tofInput.GetStrip();
374       detectorIndex[3] = (Int_t)tofInput.GetPadZ();
375       detectorIndex[4] = (Int_t)tofInput.GetPadX();
376
377       if (detectorIndex[0]==-1 ||
378           detectorIndex[1]==-1 ||
379           detectorIndex[2]==-1 ||
380           detectorIndex[3]==-1 ||
381           detectorIndex[4]==-1) continue;
382
383       //AliInfo(Form("  %2i  %1i  %2i  %1i  %2i ",detectorIndex[0],detectorIndex[1],detectorIndex[2],detectorIndex[3],detectorIndex[4]));
384
385       for (ii=0; ii<3; ii++) position[ii] =  0.;
386
387       fTOFGeometry->GetPos(detectorIndex, position);
388       
389       cylindricalPosition[0] = (Double_t)TMath::Sqrt(position[0]*position[0] + position[1]*position[1]);
390       cylindricalPosition[1] = (Double_t)TMath::ATan2(position[1], position[0]);
391       cylindricalPosition[2] = (Double_t)position[2];
392       cylindricalPosition[3] = (Double_t)tofInput.GetTofBin();
393       cylindricalPosition[4] = (Double_t)tofInput.GetToTbin();
394       tToT = tofInput.GetToTbin();
395       tTdcND = -1.;
396     
397       AliTOFcluster *tofCluster = new AliTOFcluster(cylindricalPosition, detectorIndex);
398       tofCluster->SetToT(tToT);
399       tofCluster->SetTDCND(tTdcND);
400       InsertCluster(tofCluster);
401
402     } // while loop
403
404   } // DDL Loop
405
406   AliInfo(Form("Number of found clusters: %i", fNumberOfTofClusters));
407
408   CalibrateRecPoint();
409   FillRecPoint();
410
411   fTreeR->Fill();
412   ResetRecpoint();
413
414   fTOFLoader = fRunLoader->GetLoader("TOFLoader");
415   fTOFLoader->WriteRecPoints("OVERWRITE");
416   
417 }
418 //______________________________________________________________________________
419
420 void AliTOFClusterFinder::Raw2Digits(Int_t iEvent, AliRawReader *rawReader)
421 {
422   //
423   // Converts RAW data to MC digits for TOF
424   //
425   //             (temporary solution)
426   //
427
428   fRunLoader->GetEvent(iEvent);
429
430   fTreeD = fTOFLoader->TreeD();
431   if (fTreeD)
432     {
433     AliInfo("TreeD re-creation");
434     fTreeD = 0x0;
435     fTOFLoader->MakeTree("D");
436     fTreeD = fTOFLoader->TreeD();
437     }
438
439   TClonesArray dummy("AliTOFdigit",10000), *tofDigits=&dummy;
440   Int_t bufsize = 32000;
441   fTreeD->Branch("TOF", &tofDigits, bufsize);
442
443
444   const Int_t kDDL = fTOFGeometry->NDDL()*fTOFGeometry->NSectors();
445
446   fRunLoader->GetEvent(iEvent);
447
448   AliDebug(2,Form(" Event number %2i ", iEvent));
449
450   Int_t ii = 0;
451   Int_t indexDDL = 0;
452
453   Int_t detectorIndex[5];
454   Float_t digit[4];
455
456   for (indexDDL = 0; indexDDL < kDDL; indexDDL++) {
457
458     rawReader->Reset();
459     AliTOFRawStream tofInput(rawReader);
460     rawReader->Select("TOF", indexDDL, indexDDL);
461
462     while(tofInput.Next()) {
463
464       for (ii=0; ii<5; ii++) detectorIndex[ii] = -1;
465
466       detectorIndex[0] = tofInput.GetSector();
467       detectorIndex[1] = tofInput.GetPlate();
468       detectorIndex[2] = tofInput.GetStrip();
469       detectorIndex[3] = tofInput.GetPadX();
470       detectorIndex[4] = tofInput.GetPadZ();
471       
472       //AliInfo(Form("  %2i  %1i  %2i  %1i  %2i ",detectorIndex[0],detectorIndex[1],detectorIndex[2],detectorIndex[3],detectorIndex[4]));
473
474       if (detectorIndex[0]==-1 ||
475           detectorIndex[1]==-1 ||
476           detectorIndex[2]==-1 ||
477           detectorIndex[3]==-1 ||
478           detectorIndex[4]==-1) continue;
479
480       digit[0] = (Float_t)tofInput.GetTofBin();
481       digit[1] = (Float_t)tofInput.GetToTbin();
482       digit[2] = (Float_t)tofInput.GetToTbin();
483       digit[3] = -1.;
484
485       Int_t tracknum[3]={-1,-1,-1};
486
487       TClonesArray &aDigits = *tofDigits;
488       Int_t last=tofDigits->GetEntriesFast();
489       new (aDigits[last]) AliTOFdigit(tracknum, detectorIndex, digit);
490
491     } // while loop
492
493   } // DDL Loop
494
495   fTreeD->Fill();
496
497   fTOFLoader = fRunLoader->GetLoader("TOFLoader");
498   fTOFLoader->WriteDigits("OVERWRITE");
499   
500 }
501 //______________________________________________________________________________
502
503 Int_t AliTOFClusterFinder::InsertCluster(AliTOFcluster *tofCluster) {
504   //---------------------------------------------------------------------------//
505   // This function adds a TOF cluster to the array of TOF clusters sorted in Z //
506   //---------------------------------------------------------------------------//
507   if (fNumberOfTofClusters==kTofMaxCluster) {
508     AliError("Too many clusters !");
509     return 1;
510   }
511
512   if (fNumberOfTofClusters==0) {
513     fTofClusters[fNumberOfTofClusters++] = tofCluster;
514     return 0;
515   }
516
517   Int_t ii = FindClusterIndex(tofCluster->GetZ());
518   memmove(fTofClusters+ii+1 ,fTofClusters+ii,(fNumberOfTofClusters-ii)*sizeof(AliTOFcluster*));
519   fTofClusters[ii] = tofCluster;
520   fNumberOfTofClusters++;
521   
522   return 0;
523
524 }
525 //_________________________________________________________________________
526
527 Int_t AliTOFClusterFinder::FindClusterIndex(Double_t z) const {
528   //--------------------------------------------------------------------
529   // This function returns the index of the nearest cluster 
530   //--------------------------------------------------------------------
531   if (fNumberOfTofClusters==0) return 0;
532   if (z <= fTofClusters[0]->GetZ()) return 0;
533   if (z > fTofClusters[fNumberOfTofClusters-1]->GetZ()) return fNumberOfTofClusters;
534   Int_t b = 0, e = fNumberOfTofClusters-1, m = (b+e)/2;
535   for (; b<e; m=(b+e)/2) {
536     if (z > fTofClusters[m]->GetZ()) b=m+1;
537     else e=m;
538   }
539
540   return m;
541
542 }
543 //_________________________________________________________________________
544
545 void AliTOFClusterFinder::FillRecPoint()
546 {
547   //
548   // Copy the global array of AliTOFcluster, i.e. fTofClusters (sorted
549   // in Z) in the global TClonesArray of AliTOFcluster,
550   // i.e. fRecPoints.
551   //
552
553   Int_t ii, jj;
554
555   Int_t detectorIndex[5];
556   Double_t cylindricalPosition[5];
557   Int_t trackLabels[3];
558   Int_t digitIndex = -1;
559   Float_t tToT;
560   Double_t tTdcND;
561
562   TClonesArray &lRecPoints = *fRecPoints;
563   
564   for (ii=0; ii<fNumberOfTofClusters; ii++) {
565
566     digitIndex = fTofClusters[ii]->GetIndex();
567     for(jj=0; jj<5; jj++) detectorIndex[jj] = fTofClusters[ii]->GetDetInd(jj);
568     for(jj=0; jj<3; jj++) trackLabels[jj] = fTofClusters[ii]->GetLabel(jj);
569     cylindricalPosition[0] = fTofClusters[ii]->GetR();
570     cylindricalPosition[1] = fTofClusters[ii]->GetPhi();
571     cylindricalPosition[2] = fTofClusters[ii]->GetZ();
572     cylindricalPosition[3] = fTofClusters[ii]->GetTDC();
573     cylindricalPosition[4] = fTofClusters[ii]->GetADC();
574     tToT = fTofClusters[ii]->GetToT();
575     tTdcND = fTofClusters[ii]->GetTDCND();
576
577     new(lRecPoints[ii]) AliTOFcluster(cylindricalPosition, trackLabels, detectorIndex, digitIndex, tToT, tTdcND);
578
579     //AliInfo(Form("%3i  %3i  %f %f %f %f %f  %2i %2i %2i %1i %2i",ii,digitIndex, cylindricalPosition[2],cylindricalPosition[0],cylindricalPosition[1],cylindricalPosition[3],cylindricalPosition[4],detectorIndex[0],detectorIndex[1],detectorIndex[2],detectorIndex[3],detectorIndex[4]));
580
581   } // loop on clusters
582
583 }
584
585 //_________________________________________________________________________
586 void AliTOFClusterFinder::CalibrateRecPoint()
587 {
588   //
589   // Copy the global array of AliTOFcluster, i.e. fTofClusters (sorted
590   // in Z) in the global TClonesArray of AliTOFcluster,
591   // i.e. fRecPoints.
592   //
593
594   Int_t ii, jj;
595
596   Int_t detectorIndex[5];
597   Int_t digitIndex = -1;
598   Float_t tToT;
599   Float_t tdcCorr;
600   AliInfo(" Calibrating TOF Clusters: ")
601   AliTOFcalib *calib = new AliTOFcalib(fTOFGeometry);
602   // calib->ReadParFromCDB("TOF/Calib",0); // original
603   calib->ReadParFromCDB("TOF/Calib",-1); // Use AliCDBManager's run number
604   AliTOFCal *calTOFArray = calib->GetTOFCalArray();  
605
606   for (ii=0; ii<fNumberOfTofClusters; ii++) {
607     digitIndex = fTofClusters[ii]->GetIndex();
608     for(jj=0; jj<5; jj++) detectorIndex[jj] = fTofClusters[ii]->GetDetInd(jj);
609
610     Int_t index = calib->GetIndex(detectorIndex);
611      
612     AliTOFChannel * calChannel = calTOFArray->GetChannel(index);
613     Float_t par[6];
614     for (Int_t j = 0; j<6; j++){
615       par[j]=calChannel->GetSlewPar(j);
616     }
617     tToT = fTofClusters[ii]->GetToT();
618     Float_t timeCorr=par[0]+par[1]*tToT+par[2]*tToT*tToT+par[3]*tToT*tToT*tToT+par[4]*tToT*tToT*tToT*tToT+par[5]*tToT*tToT*tToT*tToT*tToT;
619     tdcCorr=(fTofClusters[ii]->GetTDC()*AliTOFGeometry::TdcBinWidth()+32)*1.E-3-timeCorr;
620     tdcCorr=(tdcCorr*1E3-32)/AliTOFGeometry::TdcBinWidth();
621     fTofClusters[ii]->SetTDC(tdcCorr);
622
623   } // loop on clusters
624
625   delete calib;
626 }
627 //______________________________________________________________________________
628
629 void AliTOFClusterFinder::ResetRecpoint()
630 {
631   //
632   // Clear the list of reconstructed points
633   //
634
635   fNumberOfTofClusters = 0;
636   if (fRecPoints) fRecPoints->Clear();
637
638 }
639 //______________________________________________________________________________
640
641 void AliTOFClusterFinder::Load()
642 {
643   //
644   // Load TOF.Digits.root and TOF.RecPoints.root files
645   //
646
647   fTOFLoader->LoadDigits("READ");
648   fTOFLoader->LoadRecPoints("recreate");
649
650 }
651 //______________________________________________________________________________
652
653 void AliTOFClusterFinder::LoadClusters()
654 {
655   //
656   // Load TOF.RecPoints.root file
657   //
658
659   fTOFLoader->LoadRecPoints("recreate");
660
661 }
662 //______________________________________________________________________________
663
664 void AliTOFClusterFinder::UnLoad()
665 {
666   //
667   // Unload TOF.Digits.root and TOF.RecPoints.root files
668   //
669
670   fTOFLoader->UnloadDigits();
671   fTOFLoader->UnloadRecPoints();
672
673 }
674 //______________________________________________________________________________
675
676 void AliTOFClusterFinder::UnLoadClusters()
677 {
678   //
679   // Unload TOF.RecPoints.root file
680   //
681
682   fTOFLoader->UnloadRecPoints();
683
684 }
685 //______________________________________________________________________________