New detector numbering scheme (common for DAQ/HLT/Offline). All the subdetectors...
[u/mrichter/AliRoot.git] / PMD / AliPMDClusterFinder.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 //           Date   : August 05 2003                   //
19 //  This reads the file PMD.digits.root(TreeD),        //
20 //  calls the Clustering algorithm and stores the      //
21 //  clustering output in PMD.RecPoints.root(TreeR)     // 
22 //                                                     //
23 //-----------------------------------------------------//
24
25 #include <Riostream.h>
26 #include <TTree.h>
27 #include <TObjArray.h>
28 #include <TClonesArray.h>
29
30 #include "AliLog.h"
31 #include "AliRunLoader.h"
32 #include "AliLoader.h"
33 #include "AliRawReader.h"
34
35 #include "AliPMDdigit.h"
36 #include "AliPMDClusterFinder.h"
37 #include "AliPMDClustering.h"
38 #include "AliPMDClusteringV1.h"
39 #include "AliPMDcluster.h"
40 #include "AliPMDrecpoint1.h"
41 #include "AliPMDrechit.h"
42 #include "AliPMDRawStream.h"
43 #include "AliPMDCalibData.h"
44
45 #include "AliCDBManager.h"
46 #include "AliCDBEntry.h"
47
48
49
50 ClassImp(AliPMDClusterFinder)
51
52 AliPMDClusterFinder::AliPMDClusterFinder():
53   fRunLoader(0),
54   fPMDLoader(0),
55   fTreeD(0),
56   fTreeR(0),
57   fDigits(new TClonesArray("AliPMDdigit", 1000)),
58   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
59   fRechits(new TClonesArray("AliPMDrechit", 1000)),
60   fNpoint(0),
61   fNhit(0),
62   fEcut(0.)
63 {
64 //
65 // Constructor
66 //
67   fCalibData = GetCalibData();
68 }
69 // ------------------------------------------------------------------------- //
70 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
71   fRunLoader(runLoader),
72   fPMDLoader(runLoader->GetLoader("PMDLoader")),
73   fTreeD(0),
74   fTreeR(0),
75   fDigits(new TClonesArray("AliPMDdigit", 1000)),
76   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
77   fRechits(new TClonesArray("AliPMDrechit", 1000)),
78   fNpoint(0),
79   fNhit(0),
80   fEcut(0.)
81 {
82 //
83 // Constructor
84 //
85   fCalibData = GetCalibData();
86 }
87 // ------------------------------------------------------------------------- //
88 AliPMDClusterFinder::~AliPMDClusterFinder()
89 {
90   // Destructor
91   if (fDigits)
92     {
93       fDigits->Delete();
94       delete fDigits;
95       fDigits=0;
96     }
97   if (fRecpoints)
98     {
99       fRecpoints->Delete();
100       delete fRecpoints;
101       fRecpoints=0;
102     }
103   if (fRechits)
104     {
105       fRechits->Delete();
106       delete fRechits;
107       fRechits=0;
108     }
109 }
110 // ------------------------------------------------------------------------- //
111
112 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
113 {
114   // Converts digits to recpoints after running clustering
115   // algorithm on CPV plane and PREshower plane
116   //
117   Int_t    det = 0,smn = 0;
118   Int_t    xpos,ypos;
119   Float_t  adc;
120   Int_t    ismn;
121   Int_t    idet;
122   Float_t  clusdata[6];
123
124   TObjArray *pmdcont = new TObjArray();
125   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
126
127   pmdclust->SetEdepCut(fEcut);
128
129   fRunLoader->GetEvent(ievt);
130
131
132   fTreeD = fPMDLoader->TreeD();
133   if (fTreeD == 0x0)
134     {
135       AliFatal("AliPMDClusterFinder: Can not get TreeD");
136
137     }
138   AliPMDdigit  *pmddigit;
139   TBranch *branch = fTreeD->GetBranch("PMDDigit");
140   branch->SetAddress(&fDigits);
141
142   ResetRecpoint();
143
144   fTreeR = fPMDLoader->TreeR();
145   if (fTreeR == 0x0)
146     {
147       fPMDLoader->MakeTree("R");
148       fTreeR = fPMDLoader->TreeR();
149     }
150
151   Int_t bufsize = 16000;
152   TBranch * branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize); 
153   TBranch * branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize); 
154
155   Int_t nmodules = (Int_t) fTreeD->GetEntries();
156
157   for (Int_t imodule = 0; imodule < nmodules; imodule++)
158     {
159       ResetCellADC();
160       fTreeD->GetEntry(imodule); 
161       Int_t nentries = fDigits->GetLast();
162       for (Int_t ient = 0; ient < nentries+1; ient++)
163         {
164           pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
165           
166           det    = pmddigit->GetDetector();
167           smn    = pmddigit->GetSMNumber();
168           xpos   = pmddigit->GetRow();
169           ypos   = pmddigit->GetColumn();
170           adc    = pmddigit->GetADC();
171           
172           // CALIBRATION
173           Float_t gain = fCalibData->GetGainFact(det,smn,xpos,ypos);
174           
175          // printf("adc = %d gain = %f\n",adc,gain);
176           
177           adc = adc*gain;
178
179
180           //Int_t trno   = pmddigit->GetTrackNumber();
181           fCellADC[xpos][ypos] = (Double_t) adc;
182         }
183
184       idet = det;
185       ismn = smn;
186       pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
187       
188       Int_t nentries1 = pmdcont->GetEntries();
189
190       AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
191
192       for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
193         {
194           AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
195           idet        = pmdcl->GetDetector();
196           ismn        = pmdcl->GetSMN();
197           clusdata[0] = pmdcl->GetClusX();
198           clusdata[1] = pmdcl->GetClusY();
199           clusdata[2] = pmdcl->GetClusADC();
200           clusdata[3] = pmdcl->GetClusCells();
201           clusdata[4] = pmdcl->GetClusSigmaX();
202           clusdata[5] = pmdcl->GetClusSigmaY();
203
204           AddRecPoint(idet,ismn,clusdata);
205
206           Int_t ncell = (Int_t) clusdata[3];
207           for(Int_t ihit = 0; ihit < ncell; ihit++)
208             {
209               Int_t celldataX = pmdcl->GetClusCellX(ihit);
210               Int_t celldataY = pmdcl->GetClusCellY(ihit);
211               AddRecHit(celldataX, celldataY);
212             }
213           branch2->Fill();
214           ResetRechit();
215         }
216       pmdcont->Clear();
217       
218       branch1->Fill();
219       ResetRecpoint();
220
221     } // modules
222
223   ResetCellADC();
224   fPMDLoader = fRunLoader->GetLoader("PMDLoader");  
225   fPMDLoader->WriteRecPoints("OVERWRITE");
226
227   //   delete the pointers
228   delete pmdclust;
229   delete pmdcont;
230     
231 }
232 // ------------------------------------------------------------------------- //
233
234 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
235                                            TTree *clustersTree)
236 {
237   // Converts RAW data to recpoints after running clustering
238   // algorithm on CPV and PREshower plane
239   //
240
241   Float_t  clusdata[6];
242
243   TObjArray *pmdcont = new TObjArray();
244   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
245
246   pmdclust->SetEdepCut(fEcut);
247
248   ResetRecpoint();
249
250   Int_t bufsize = 16000;
251   TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize); 
252
253   TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize); 
254
255   const Int_t kDDL = 6;
256   const Int_t kRow = 48;
257   const Int_t kCol = 96;
258
259   Int_t idet = 0;
260   Int_t iSMN = 0;
261   
262   for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
263     {
264       if (indexDDL < 4)
265         {
266           iSMN = 6;
267         }
268       else if (indexDDL >= 4)
269         {
270           iSMN = 12;
271         }
272       Int_t ***precpvADC;
273       precpvADC = new int **[iSMN];
274       for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
275       for (Int_t i=0; i<iSMN;i++)
276         {
277           for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
278         }
279       for (Int_t i = 0; i < iSMN; i++)
280         {
281           for (Int_t j = 0; j < kRow; j++)
282             {
283               for (Int_t k = 0; k < kCol; k++)
284                 {
285                   precpvADC[i][j][k] = 0;
286                 }
287             }
288         }
289       ResetCellADC();
290       rawReader->Reset();
291       AliPMDRawStream pmdinput(rawReader);
292       rawReader->Select("PMD", indexDDL, indexDDL);
293       while(pmdinput.Next())
294         {
295           Int_t det = pmdinput.GetDetector();
296           Int_t smn = pmdinput.GetSMN();
297           //Int_t mcm = pmdinput.GetMCM();
298           //Int_t chno = pmdinput.GetChannel();
299           Int_t row = pmdinput.GetRow();
300           Int_t col = pmdinput.GetColumn();
301           Int_t sig = pmdinput.GetSignal();
302           
303           Float_t sig1 = (Float_t) sig;
304
305           // CALIBRATION
306           Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
307           
308         //  printf("sig = %d gain = %f\n",sig,gain);
309           
310           sig = (Int_t) (sig1*gain);
311           
312           Int_t indexsmn = 0;
313
314           if (indexDDL < 4)
315             {
316               if (det != 0)
317                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
318                               indexDDL, det));
319               indexsmn = smn - indexDDL * 6;
320             }
321           else if (indexDDL == 4)
322             {
323               if (det != 1)
324                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
325                               indexDDL, det));
326               if (smn < 6)
327                 {
328                   indexsmn = smn;
329                 }
330               else if (smn >= 12 && smn < 18)
331                 {
332                   indexsmn = smn - 6;
333                 }
334             }
335           else if (indexDDL == 5)
336             {
337               if (det != 1)
338                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
339                               indexDDL, det));
340               if (smn >= 6 && smn < 12)
341                 {
342                   indexsmn = smn - 6;
343                 }
344               else if (smn >= 18 && smn < 24)
345                 {
346                   indexsmn = smn - 12;
347                 }
348             }         
349           precpvADC[indexsmn][row][col] = sig;
350         } // while loop
351
352       Int_t ismn = 0;
353       for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
354         {
355           ResetCellADC();
356           for (Int_t irow = 0; irow < kRow; irow++)
357             {
358               for (Int_t icol = 0; icol < kCol; icol++)
359                 {
360                   fCellADC[irow][icol] = 
361                     (Double_t) precpvADC[indexsmn][irow][icol];
362                 } // row
363             }     // col
364           if (indexDDL < 4)
365             {
366               ismn = indexsmn + indexDDL * 6;
367               idet = 0;
368             }
369           else if (indexDDL == 4)
370             {
371               if (indexsmn < 6)
372                 {
373                   ismn = indexsmn;
374                 }
375               else if (indexsmn >= 6 && indexsmn < 12)
376                 {
377                   ismn = indexsmn + 6;
378                 }
379               idet = 1;
380             }
381           else if (indexDDL == 5)
382             {
383               if (indexsmn < 6)
384                 {
385                   ismn = indexsmn + 6;
386                 }
387               else if (indexsmn >= 6 && indexsmn < 12)
388                 {
389                   ismn = indexsmn + 12;
390                 }
391               idet = 1;
392             }
393
394
395           pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
396           Int_t nentries1 = pmdcont->GetEntries();
397
398           AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
399
400           for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
401             {
402               AliPMDcluster *pmdcl = 
403                 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
404               idet        = pmdcl->GetDetector();
405               ismn        = pmdcl->GetSMN();
406               clusdata[0] = pmdcl->GetClusX();
407               clusdata[1] = pmdcl->GetClusY();
408               clusdata[2] = pmdcl->GetClusADC();
409               clusdata[3] = pmdcl->GetClusCells();
410               clusdata[4] = pmdcl->GetClusSigmaX();
411               clusdata[5] = pmdcl->GetClusSigmaY();
412
413               AddRecPoint(idet,ismn,clusdata);
414
415               Int_t ncell = (Int_t) clusdata[3];
416               for(Int_t ihit = 0; ihit < ncell; ihit++)
417                 {
418                   Int_t celldataX = pmdcl->GetClusCellX(ihit);
419                   Int_t celldataY = pmdcl->GetClusCellY(ihit);
420                   AddRecHit(celldataX, celldataY);
421                 }
422               branch2->Fill();
423               ResetRechit();
424
425             }
426           pmdcont->Clear();
427           
428           branch1->Fill();
429           ResetRecpoint();
430
431
432         } // smn
433
434       for (Int_t i=0; i<iSMN; i++)
435         {
436           for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
437         }
438       for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
439       delete precpvADC;
440     } // DDL Loop
441   
442   ResetCellADC();
443   
444   //   delete the pointers
445   delete pmdclust;
446   delete pmdcont;
447
448 }
449 // ------------------------------------------------------------------------- //
450
451 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
452 {
453   // Converts RAW data to recpoints after running clustering
454   // algorithm on CPV and PREshower plane
455   //
456
457   Float_t  clusdata[6];
458
459   TObjArray *pmdcont = new TObjArray();
460
461   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
462
463   pmdclust->SetEdepCut(fEcut);
464
465   fRunLoader->GetEvent(ievt);
466
467   ResetRecpoint();
468
469   fTreeR = fPMDLoader->TreeR();
470   if (fTreeR == 0x0)
471     {
472       fPMDLoader->MakeTree("R");
473       fTreeR = fPMDLoader->TreeR();
474     }
475   Int_t bufsize = 16000;
476   TBranch *branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize); 
477   TBranch *branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize); 
478
479   const Int_t kDDL = 6;
480   const Int_t kRow = 48;
481   const Int_t kCol = 96;
482
483   Int_t idet = 0;
484   Int_t iSMN = 0;
485   
486   for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
487     {
488       if (indexDDL < 4)
489         {
490           iSMN = 6;
491         }
492       else if (indexDDL >= 4)
493         {
494           iSMN = 12;
495         }
496       Int_t ***precpvADC;
497       precpvADC = new int **[iSMN];
498       for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
499       for (Int_t i=0; i<iSMN;i++)
500         {
501           for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
502         }
503       for (Int_t i = 0; i < iSMN; i++)
504         {
505           for (Int_t j = 0; j < kRow; j++)
506             {
507               for (Int_t k = 0; k < kCol; k++)
508                 {
509                   precpvADC[i][j][k] = 0;
510                 }
511             }
512         }
513       ResetCellADC();
514       rawReader->Reset();
515       AliPMDRawStream pmdinput(rawReader);
516       rawReader->Select("PMD", indexDDL, indexDDL);
517       while(pmdinput.Next())
518         {
519           Int_t det = pmdinput.GetDetector();
520           Int_t smn = pmdinput.GetSMN();
521           //Int_t mcm = pmdinput.GetMCM();
522           //Int_t chno = pmdinput.GetChannel();
523           Int_t row = pmdinput.GetRow();
524           Int_t col = pmdinput.GetColumn();
525           Int_t sig = pmdinput.GetSignal();
526
527           Float_t sig1 = (Float_t) sig;
528           // CALIBRATION
529           Float_t gain = fCalibData->GetGainFact(det,smn,row,col);
530           
531          // printf("sig = %d gain = %f\n",sig,gain);
532           
533           sig = (Int_t) (sig1*gain);
534           
535           Int_t indexsmn = 0;
536
537           if (indexDDL < 4)
538             {
539               if (det != 0)
540                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
541                               indexDDL, det));
542               indexsmn = smn - indexDDL * 6;
543             }
544           else if (indexDDL == 4)
545             {
546               if (det != 1)
547                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
548                               indexDDL, det));
549               if (smn < 6)
550                 {
551                   indexsmn = smn;
552                 }
553               else if (smn >= 12 && smn < 18)
554                 {
555                   indexsmn = smn - 6;
556                 }
557             }
558           else if (indexDDL == 5)
559             {
560               if (det != 1)
561                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
562                               indexDDL, det));
563               if (smn >= 6 && smn < 12)
564                 {
565                   indexsmn = smn - 6;
566                 }
567               else if (smn >= 18 && smn < 24)
568                 {
569                   indexsmn = smn - 12;
570                 }
571             }         
572           precpvADC[indexsmn][row][col] = sig;
573         } // while loop
574
575       Int_t ismn = 0;
576       for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
577         {
578           ResetCellADC();
579           for (Int_t irow = 0; irow < kRow; irow++)
580             {
581               for (Int_t icol = 0; icol < kCol; icol++)
582                 {
583                   fCellADC[irow][icol] = 
584                     (Double_t) precpvADC[indexsmn][irow][icol];
585                 } // row
586             }     // col
587           if (indexDDL < 4)
588             {
589               ismn = indexsmn + indexDDL * 6;
590               idet = 0;
591             }
592           else if (indexDDL == 4)
593             {
594               if (indexsmn < 6)
595                 {
596                   ismn = indexsmn;
597                 }
598               else if (indexsmn >= 6 && indexsmn < 12)
599                 {
600                   ismn = indexsmn + 6;
601                 }
602               idet = 1;
603             }
604           else if (indexDDL == 5)
605             {
606               if (indexsmn < 6)
607                 {
608                   ismn = indexsmn + 6;
609                 }
610               else if (indexsmn >= 6 && indexsmn < 12)
611                 {
612                   ismn = indexsmn + 12;
613                 }
614               idet = 1;
615             }
616
617
618           pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
619           Int_t nentries1 = pmdcont->GetEntries();
620
621           AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
622
623           for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
624             {
625               AliPMDcluster *pmdcl = 
626                 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
627               idet        = pmdcl->GetDetector();
628               ismn        = pmdcl->GetSMN();
629               clusdata[0] = pmdcl->GetClusX();
630               clusdata[1] = pmdcl->GetClusY();
631               clusdata[2] = pmdcl->GetClusADC();
632               clusdata[3] = pmdcl->GetClusCells();
633               clusdata[4] = pmdcl->GetClusSigmaX();
634               clusdata[5] = pmdcl->GetClusSigmaY();
635
636               AddRecPoint(idet,ismn,clusdata);
637
638               Int_t ncell = (Int_t) clusdata[3];
639               for(Int_t ihit = 0; ihit < ncell; ihit++)
640                 {
641                   Int_t celldataX = pmdcl->GetClusCellX(ihit);
642                   Int_t celldataY = pmdcl->GetClusCellY(ihit);
643                   AddRecHit(celldataX, celldataY);
644                 }
645               branch2->Fill();
646               ResetRechit();
647
648             }
649           pmdcont->Clear();
650           
651           branch1->Fill();
652           ResetRecpoint();
653
654
655         } // smn
656
657       for (Int_t i=0; i<iSMN; i++)
658         {
659           for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
660         }
661       for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
662       delete precpvADC;
663     } // DDL Loop
664   
665   ResetCellADC();
666   
667   fPMDLoader = fRunLoader->GetLoader("PMDLoader");  
668   fPMDLoader->WriteRecPoints("OVERWRITE");
669
670   //   delete the pointers
671   delete pmdclust;
672   delete pmdcont;
673
674 }
675 // ------------------------------------------------------------------------- //
676 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
677 {
678   fEcut = ecut;
679 }
680 // ------------------------------------------------------------------------- //
681 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
682 {
683   // Add Reconstructed points
684   //
685   TClonesArray &lrecpoints = *fRecpoints;
686   AliPMDrecpoint1 *newrecpoint;
687   newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
688   new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
689   delete newrecpoint;
690 }
691 // ------------------------------------------------------------------------- //
692 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY)
693 {
694   // Add associated cell hits to the Reconstructed points
695   //
696   TClonesArray &lrechits = *fRechits;
697   AliPMDrechit *newrechit;
698   newrechit = new AliPMDrechit(celldataX, celldataY);
699   new(lrechits[fNhit++]) AliPMDrechit(newrechit);
700   delete newrechit;
701 }
702 // ------------------------------------------------------------------------- //
703 void AliPMDClusterFinder::ResetCellADC()
704 {
705   // Reset the individual cell ADC value to zero
706   //
707   for(Int_t irow = 0; irow < fgkRow; irow++)
708     {
709       for(Int_t icol = 0; icol < fgkCol; icol++)
710         {
711           fCellADC[irow][icol] = 0.;
712         }
713     }
714 }
715 // ------------------------------------------------------------------------- //
716
717 void AliPMDClusterFinder::ResetRecpoint()
718 {
719   // Clear the list of reconstructed points
720   fNpoint = 0;
721   if (fRecpoints) fRecpoints->Clear();
722 }
723 // ------------------------------------------------------------------------- //
724 void AliPMDClusterFinder::ResetRechit()
725 {
726   // Clear the list of reconstructed points
727   fNhit = 0;
728   if (fRechits) fRechits->Clear();
729 }
730 // ------------------------------------------------------------------------- //
731 void AliPMDClusterFinder::Load()
732 {
733   // Load all the *.root files
734   //
735   fPMDLoader->LoadDigits("READ");
736   fPMDLoader->LoadRecPoints("recreate");
737 }
738 // ------------------------------------------------------------------------- //
739 void AliPMDClusterFinder::LoadClusters()
740 {
741   // Load all the *.root files
742   //
743   fPMDLoader->LoadRecPoints("recreate");
744 }
745 // ------------------------------------------------------------------------- //
746 void AliPMDClusterFinder::UnLoad()
747 {
748   // Unload all the *.root files
749   //
750   fPMDLoader->UnloadDigits();
751   fPMDLoader->UnloadRecPoints();
752 }
753 // ------------------------------------------------------------------------- //
754 void AliPMDClusterFinder::UnLoadClusters()
755 {
756   // Unload all the *.root files
757   //
758   fPMDLoader->UnloadRecPoints();
759 }
760 // ------------------------------------------------------------------------- //
761
762 AliPMDCalibData* AliPMDClusterFinder::GetCalibData() const
763 {
764   // The run number will be centralized in AliCDBManager,
765   // you don't need to set it here!
766   // Added by ZA
767   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Data");
768   
769   if(!entry){
770     AliWarning("Calibration object retrieval failed! Dummy calibration will be used.");
771     
772     // this just remembers the actual default storage. No problem if it is null.
773     AliCDBStorage *origStorage = AliCDBManager::Instance()->GetDefaultStorage();
774     AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
775     
776     entry = AliCDBManager::Instance()->Get("PMD/Calib/Data");
777     
778     // now reset the original default storage to AliCDBManager...
779     AliCDBManager::Instance()->SetDefaultStorage(origStorage);  
780   }
781   
782   AliPMDCalibData *calibdata=0;
783   if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
784   
785   if (!calibdata)  AliError("No calibration data from calibration database !");
786   
787   return calibdata;
788 }