1899a5ab0aa75c0f1493f86d3f6718db1b70ccf8
[u/mrichter/AliRoot.git] / HLT / src / AliL3FileHandler.cxx
1 //$Id$
2
3 // Author: Uli Frankenfeld <mailto:franken@fi.uib.no>, Anders Vestbo <mailto:vestbo$fi.uib.no>
4 //*-- Copyright &copy Uli 
5
6 #include "AliL3StandardIncludes.h"
7 #include <AliTPCDigitsArray.h>
8 #include <AliTPCClustersArray.h>
9 #include <AliTPCcluster.h>
10 #include <AliTPCClustersRow.h>
11
12 #include "AliL3Logging.h"
13 #include "AliL3Transform.h"
14 #include "AliL3MemHandler.h"
15 #include "AliL3FileHandler.h"
16
17 #include "AliL3DigitData.h"
18 #include "AliL3TrackSegmentData.h"
19 #include "AliL3SpacePointData.h"
20 #include "AliL3TrackArray.h"
21
22 //_____________________________________________________________
23 // AliL3FileHandler
24 //
25 // The HLT ROOT <-> binary files handling class
26 //
27 // This class provides the interface between AliROOT files,
28 // and HLT binary files. It should be used for converting 
29 // TPC data stored in AliROOT format (outputfile from a simulation),
30 // into the data format currently used by in the HLT framework. 
31 // This enables the possibility to always use the same data format, 
32 // whether you are using a binary file as an input, or a AliROOT file.
33 //
34 // For example on how to create binary files from a AliROOT simulation,
35 // see example macro exa/Binary.C.
36 //
37 // For reading a AliROOT file into HLT format in memory, do the following:
38 //
39 // AliL3FileHandler file;
40 // file.SetAliInput("galice.root");
41 // AliL3DigitRowData *dataPt = (AliL3DigitRowData*)file.AliDigits2Memory(nrows,eventnr);
42 // 
43 // All the data are then stored in memory and accessible via the pointer dataPt.
44 // Accesing the data is then identical to the example 1) showed in AliL3MemHandler class.
45 //
46 // For converting the data back, and writing it to a new AliROOT file do:
47 //
48 // AliL3FileHandler file;
49 // file.SetAliInput("galice.root");
50 // file.Init(slice,patch,NumberOfRowsInPatch);
51 // file.AliDigits2RootFile(dataPt,"new_galice.root");
52 // file.CloseAliInput();
53
54 ClassImp(AliL3FileHandler)
55
56 AliL3FileHandler::AliL3FileHandler()
57 {
58   //Default constructor
59   fInAli = 0;
60   fParam = 0;
61   fMC =0;
62   fLastIndex=0;
63   fDigits=0;
64   fDigitsTree=0;
65 }
66
67 AliL3FileHandler::~AliL3FileHandler()
68 {
69   //Destructor
70   if(fMC) CloseMCOutput();
71   FreeDigitsTree();
72   if(fInAli) CloseAliInput();
73   
74 }
75
76 void AliL3FileHandler::FreeDigitsTree()
77 {
78   if(!fDigitsTree)
79     {
80       LOG(AliL3Log::kWarning,"AliL3FileHandler::FreeDigitsTree()","Pointer")
81         <<"Cannot free digitstree, it is not present"<<ENDLOG;
82       return;
83     }
84   fDigits=0;
85   fDigitsTree->Delete();
86   fDigitsTree=0;
87 }
88
89 Bool_t AliL3FileHandler::SetMCOutput(char *name)
90 {
91   fMC = fopen(name,"w");
92   if(!fMC){
93     LOG(AliL3Log::kWarning,"AliL3FileHandler::SetMCOutput","File Open")
94       <<"Pointer to File = 0x0 "<<ENDLOG;
95     return kFALSE;
96   }
97   return kTRUE;
98 }
99
100 Bool_t AliL3FileHandler::SetMCOutput(FILE *file)
101 {
102   fMC = file;
103   if(!fMC){
104     LOG(AliL3Log::kWarning,"AliL3FileHandler::SetMCOutput","File Open")
105       <<"Pointer to File = 0x0 "<<ENDLOG;
106     return kFALSE;
107   }
108   return kTRUE;
109 }
110
111 void AliL3FileHandler::CloseMCOutput()
112 {
113   if(!fMC){
114     LOG(AliL3Log::kWarning,"AliL3FileHandler::CloseMCOutPut","File Close")
115       <<"Nothing to Close"<<ENDLOG;
116     return;
117   }
118   fclose(fMC);
119   fMC =0;
120 }
121
122 Bool_t AliL3FileHandler::SetAliInput()
123 {
124   if(!fInAli->IsOpen()){
125     LOG(AliL3Log::kError,"AliL3FileHandler::SetAliInput","File Open")
126       <<"Ali File "<<fInAli->GetName()<<" does not exist"<<ENDLOG;
127     return kFALSE;
128   }
129   fParam = (AliTPCParam*)fInAli->Get("75x40_100x60_150x60");
130   if(!fParam){ 
131     LOG(AliL3Log::kError,"AliL3FileHandler::SetAliInput","File Open")
132       <<"No AliTPCParam 75x40_100x60 in File "<<fInAli->GetName()<<ENDLOG;
133     return kFALSE;
134   }
135   return kTRUE;
136 }
137
138 Bool_t AliL3FileHandler::SetAliInput(char *name)
139 {
140   //Open the AliROOT file with name.
141   
142   fInAli= new TFile(name,"READ");
143   if(!fInAli){
144     LOG(AliL3Log::kWarning,"AliL3FileHandler::SetAliInput","File Open")
145     <<"Pointer to TFile = 0x0 "<<ENDLOG;
146     return kFALSE;
147   }
148   return SetAliInput();
149 }
150
151 Bool_t AliL3FileHandler::SetAliInput(TFile *file)
152 {
153   //Specify already opened AliROOT file to use as an input.
154   
155   fInAli=file;
156   if(!fInAli){
157     LOG(AliL3Log::kWarning,"AliL3FileHandler::SetAliInput","File Open")
158     <<"Pointer to TFile = 0x0 "<<ENDLOG;
159     return kFALSE;
160   }
161   return SetAliInput();
162 }
163
164 void AliL3FileHandler::CloseAliInput()
165 {
166   if(!fInAli){
167     LOG(AliL3Log::kWarning,"AliL3FileHandler::CloseAliInput","File Close")
168       <<"Nothing to Close"<<ENDLOG;
169     return;
170   }
171   if(fInAli->IsOpen()) fInAli->Close();
172   delete fInAli;
173   fInAli = 0;
174   
175 }
176
177 Bool_t AliL3FileHandler::IsDigit(Int_t event)
178 {
179   //Check if there is a TPC digit tree in the current file.
180   //Return kTRUE if tree was found, and kFALSE if not found.
181   
182   if(!fInAli){
183     LOG(AliL3Log::kWarning,"AliL3FileHandler::IsDigit","File")
184     <<"Pointer to TFile = 0x0 "<<ENDLOG;
185     return kTRUE;  //may you are use binary input which is Digits!!
186   }
187   Char_t name[1024];
188   sprintf(name,"TreeD_75x40_100x60_150x60_%d",event);
189   TTree *t=(TTree*)fInAli->Get(name);
190   if(t){
191     LOG(AliL3Log::kInformational,"AliL3FileHandler::IsDigit","File Type")
192     <<"Found Digit Tree -> Use Fast Cluster Finder"<<ENDLOG;
193     return kTRUE;
194   }
195   else{
196     LOG(AliL3Log::kInformational,"AliL3FileHandler::IsDigit","File Type")
197     <<"No Digit Tree -> Use Cluster Tree"<<ENDLOG;
198     return kFALSE;
199   }
200 }
201
202 ///////////////////////////////////////// Digit IO  
203 Bool_t AliL3FileHandler::AliDigits2Binary(Int_t event)
204 {
205   
206   Bool_t out = kTRUE;
207   UInt_t nrow;
208   AliL3DigitRowData* data = AliDigits2Memory(nrow,event);
209   out = Memory2Binary(nrow,data);
210   Free();
211   return out;
212 }
213
214 Bool_t AliL3FileHandler::AliDigits2CompBinary(Int_t event)
215 {
216   //Convert AliROOT TPC data, into HLT data format.
217   //event specifies the event you want in the aliroot file.
218   
219   Bool_t out = kTRUE;
220   UInt_t ndigits=0;
221   AliL3DigitRowData* digits = AliDigits2Memory(ndigits,event);
222   out = Memory2CompBinary(ndigits,digits);
223   Free();
224   return out;
225 }
226
227 AliL3DigitRowData * AliL3FileHandler::AliDigits2Memory(UInt_t & nrow,Int_t event)
228 {
229   //Read data from AliROOT file into memory, and store it in the HLT data format.
230   //Returns a pointer to the data.
231
232   AliL3DigitRowData *data = 0;
233   nrow=0;
234   
235   if(!fInAli){
236     LOG(AliL3Log::kWarning,"AliL3FileHandler::AliDigits2Memory","File")
237     <<"No Input avalible: no object TFile"<<ENDLOG;
238     return 0; 
239   }
240   if(!fInAli->IsOpen()){
241     LOG(AliL3Log::kWarning,"AliL3FileHandler::AliDigits2Memory","File")
242     <<"No Input avalible: TFile not opend"<<ENDLOG;
243     return 0;
244   }
245
246   if(!fDigitsTree)
247     GetDigitsTree(event);
248   
249   UShort_t dig;
250   Int_t time,pad,sector,row;
251   Int_t nrows=0;
252   Int_t ndigitcount=0;
253   Int_t entries = (Int_t)fDigitsTree->GetEntries();
254   Int_t ndigits[entries];
255   Int_t lslice,lrow;
256   Float_t xyz[3];
257   
258   for(Int_t n=fLastIndex; n<fDigitsTree->GetEntries(); n++)
259     {
260       fDigitsTree->GetEvent(n);
261       fParam->AdjustSectorRow(fDigits->GetID(),sector,row);
262       AliL3Transform::Sector2Slice(lslice,lrow,sector,row);
263       //if(fSlice != lslice || lrow<fRowMin || lrow>fRowMax) continue;
264       if(lslice < fSlice) continue;
265       if(lslice != fSlice) break;
266       if(lrow < fRowMin) continue;
267       if(lrow > fRowMax) break;
268
269       ndigits[lrow] = 0;
270       fDigits->First();
271       do {
272         time=fDigits->CurrentRow();
273         pad=fDigits->CurrentColumn();
274         dig = fDigits->GetDigit(time,pad);
275         if(dig<=fParam->GetZeroSup()) continue;
276         
277         AliL3Transform::Raw2Local(xyz,sector,row,pad,time);
278         if(fParam->GetPadRowRadii(sector,row)<230./250.*fabs(xyz[2]))
279           continue;
280
281         ndigits[lrow]++; //for this row only
282         ndigitcount++;   //total number of digits to be published
283
284       } while (fDigits->Next());
285       //cout << lrow << " " << ndigits[lrow] << " - " << ndigitcount << endl;
286       nrows++;
287     }
288
289   Int_t size = sizeof(AliL3DigitData)*ndigitcount
290     + nrows*sizeof(AliL3DigitRowData);
291
292   LOG(AliL3Log::kDebug,"AliL3FileHandler::AliDigits2Memory","Digits")
293     <<AliL3Log::kDec<<"Found "<<ndigitcount<<" Digits"<<ENDLOG;
294   
295   data=(AliL3DigitRowData*) Allocate(size);
296   nrow = (UInt_t)nrows;
297   AliL3DigitRowData *tempPt = data;
298   for(Int_t n=fLastIndex; n<fDigitsTree->GetEntries(); n++)
299     {
300       fDigitsTree->GetEvent(n);
301       fParam->AdjustSectorRow(fDigits->GetID(),sector,row);
302       AliL3Transform::Sector2Slice(lslice,lrow,sector,row);
303       //if(fSlice != lslice || lrow<fRowMin || lrow>fRowMax) continue;
304       if(lslice < fSlice) continue;
305       if(lslice != fSlice) break;
306       if(lrow < fRowMin) continue;
307       if(lrow > fRowMax) break;
308
309       tempPt->fRow = lrow;
310       tempPt->fNDigit = ndigits[lrow];
311
312       Int_t localcount=0;
313       fDigits->First();
314       do {
315         time=fDigits->CurrentRow();
316         pad=fDigits->CurrentColumn();
317         dig = fDigits->GetDigit(time,pad);
318         if (dig <= fParam->GetZeroSup()) continue;
319         
320         //Exclude data outside cone:
321         AliL3Transform::Raw2Local(xyz,sector,row,pad,time);
322         if(fParam->GetPadRowRadii(sector,row)<230./250.*fabs(xyz[2]))
323           continue;
324
325         if(localcount >= ndigits[lrow])
326           LOG(AliL3Log::kFatal,"AliL3FileHandler::AliDigits2Binary","Memory")
327             <<AliL3Log::kDec<<"Mismatch: localcount "<<localcount<<" ndigits "
328             <<ndigits[lrow]<<ENDLOG;
329         
330         tempPt->fDigitData[localcount].fCharge=dig;
331         tempPt->fDigitData[localcount].fPad=pad;
332         tempPt->fDigitData[localcount].fTime=time;
333 #ifdef do_mc
334         tempPt->fDigitData[localcount].fTrackID[0] = fDigits->GetTrackID(time,pad,0);
335         tempPt->fDigitData[localcount].fTrackID[1] = fDigits->GetTrackID(time,pad,1);
336         tempPt->fDigitData[localcount].fTrackID[2] = fDigits->GetTrackID(time,pad,2);
337 #endif
338         localcount++;
339       } while (fDigits->Next());
340
341       Byte_t *tmp = (Byte_t*)tempPt;
342       Int_t size = sizeof(AliL3DigitRowData)
343                                       + ndigits[lrow]*sizeof(AliL3DigitData);
344       tmp += size;
345       tempPt = (AliL3DigitRowData*)tmp;
346       //fLastIndex=n;
347     }
348   //fLastIndex++;
349   return data;
350 }
351
352 Bool_t AliL3FileHandler::GetDigitsTree(Int_t event)
353 {
354   //Connects to the TPC digit tree in the AliROOT file.
355   
356   fInAli->cd();
357   Char_t dname[100];
358   sprintf(dname,"TreeD_75x40_100x60_150x60_%d",event);
359   fDigitsTree = (TTree*)fInAli->Get(dname);
360   if(!fDigitsTree) 
361     {
362       LOG(AliL3Log::kError,"AliL3FileHandler::GetDigitsTree","Digits Tree")
363         <<AliL3Log::kHex<<"Error getting digitstree "<<(Int_t)fDigitsTree<<ENDLOG;
364       return kFALSE;
365     }
366   fDigitsTree->GetBranch("Segment")->SetAddress(&fDigits);
367   return kTRUE;
368 }
369
370 void AliL3FileHandler::AliDigits2RootFile(AliL3DigitRowData *rowPt,Char_t *new_digitsfile)
371 {
372   //Write the data stored in rowPt, into a new AliROOT file.
373   //The data is stored in the AliROOT format 
374   //This is specially a nice thing if you have modified data, and wants to run it  
375   //through the offline reconstruction chain.
376   //The arguments is a pointer to the data, and the name of the new AliROOT file.
377   //Remember to pass the original AliROOT file (the one that contains the original
378   //simulated data) to this object, in order to retrieve the MC id's of the digits.
379   
380   if(!fInAli)
381     {
382       printf("AliL3FileHandler::AliDigits2RootFile : No rootfile\n");
383       return;
384     }
385   if(!fParam)
386     {
387       printf("AliL3FileHandler::AliDigits2RootFile : No parameter object. Run on rootfile\n");
388       return;
389     }
390   
391   //Get the original digitstree:
392   fInAli->cd();
393   AliTPCDigitsArray *old_array = new AliTPCDigitsArray();
394   old_array->Setup(fParam);
395   old_array->SetClass("AliSimDigits");
396   Bool_t ok = old_array->ConnectTree("TreeD_75x40_100x60_0");
397   if(!ok)
398     {
399       printf("AliL3FileHandler::AliDigits2RootFile : No digits tree object\n");
400       return;
401     }
402   
403   Bool_t create=kFALSE;
404   TFile *digFile;
405   
406   digFile = TFile::Open(new_digitsfile,"NEW");
407   if(digFile->IsOpen())
408     {    
409       create = kTRUE;
410       fParam->Write(fParam->GetTitle());
411     }
412   else
413     {
414       LOG(AliL3Log::kDebug,"AliL3FileHandler::AliDigits2RootFile","Rootfile")
415         <<"Rootfile did already exist, so I will just open it for updates"<<ENDLOG;
416       digFile = TFile::Open(new_digitsfile,"UPDATE");
417       create=kFALSE;
418     }
419   if(!digFile->IsOpen())
420     {
421       LOG(AliL3Log::kError,"AliL3FileHandler::AliDigits2RootFile","Rootfile")
422         <<"Error opening rootfile "<<new_digitsfile<<ENDLOG;
423       return;
424     }
425   
426   digFile->cd();
427     
428   //setup a new one, or connect it to the existing one:
429   AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
430   arr->SetClass("AliSimDigits");
431   arr->Setup(fParam);
432   if(create)
433     arr->MakeTree();
434   else
435     {
436       Bool_t ok = arr->ConnectTree("TreeD_75x40_100x60_0");
437       if(!ok)
438         {
439           printf("AliL3FileHandler::AliDigits2RootFile : No digits tree object in existing file\n");
440           return;
441         }
442     }
443   Int_t digcounter=0;
444
445   for(Int_t i=fRowMin; i<=fRowMax; i++)
446     {
447       
448       if((Int_t)rowPt->fRow != i) printf("AliL3FileHandler::AliDigits2RootFile : Mismatching row numbering!!!\n");
449             
450       Int_t sector,row;
451       AliL3Transform::Slice2Sector(fSlice,i,sector,row);
452       AliSimDigits * dig = (AliSimDigits*)arr->CreateRow(sector,row);
453       AliSimDigits *old_dig = (AliSimDigits*)old_array->LoadRow(sector,row);
454       if(!old_dig)
455         printf("AliL3FileHandler::AliDigits2RootFile : No padrow %d %d\n",sector,row);
456       
457       AliL3DigitData *digPt = rowPt->fDigitData;
458       digcounter=0;
459       for(UInt_t j=0; j<rowPt->fNDigit; j++)
460         {
461           Int_t charge = (Int_t)digPt[j].fCharge;
462           Int_t pad = (Int_t)digPt[j].fPad;
463           Int_t time = (Int_t)digPt[j].fTime;
464           
465           if(charge == 0) //Only write the digits that has not been removed
466             continue;
467           digcounter++;
468           dig->SetDigitFast(charge,time,pad);
469           
470           Int_t trackID[3] = {old_dig->GetTrackID(time,pad,0),old_dig->GetTrackID(time,pad,1),old_dig->GetTrackID(time,pad,2)};
471           Int_t s_pad = pad;
472           Int_t s_time = time - 1;
473           while(trackID[0] < 0)
474             {
475               if(s_time >= 0 && s_time < AliL3Transform::GetNTimeBins() && s_pad >= 0 && s_pad < AliL3Transform::GetNPads(i))
476                 {
477                   if(old_dig->GetTrackID(s_time,s_pad,0) > 0)
478                     {
479                       trackID[0]=old_dig->GetTrackID(s_time,s_pad,0); 
480                       trackID[1]=old_dig->GetTrackID(s_time,s_pad,1); 
481                       trackID[2]=old_dig->GetTrackID(s_time,s_pad,2); 
482                     }
483                 }
484               if(s_pad == pad && s_time == time - 1)
485                 s_time = time + 1;
486               else if(s_pad == pad && s_time == time + 1)
487                 {s_pad = pad - 1; s_time = time;}
488               else if(s_pad == pad - 1 && s_time == time)
489                 s_time = time - 1;
490               else if(s_pad == pad - 1 && s_time == time - 1)
491                 s_time = time + 1;
492               else if(s_pad == pad - 1 && s_time == time + 1)
493                 {s_pad = pad + 1; s_time = time;}
494               else if(s_pad == pad + 1 && s_time == time)
495                 s_time = time - 1;
496               else if(s_pad == pad + 1 && s_time == time - 1)
497                 s_time = time + 1;
498               else 
499                 break;
500             }
501           
502           dig->SetTrackIDFast(trackID[0],time,pad,0);
503           dig->SetTrackIDFast(trackID[1],time,pad,1);
504           dig->SetTrackIDFast(trackID[2],time,pad,2);
505           
506         }
507       //cout<<"Wrote "<<digcounter<<" on row "<<i<<endl;
508       UpdateRowPointer(rowPt);
509       arr->StoreRow(sector,row);
510       arr->ClearRow(sector,row);  
511       old_array->ClearRow(sector,row);
512     }
513   digFile->cd();
514   char treeName[100];
515   sprintf(treeName,"TreeD_%s_0",fParam->GetTitle());
516   printf("Writing tree to file.....");
517   arr->GetTree()->Write(treeName,TObject::kOverwrite);
518   printf("done\n");
519   digFile->Close();
520   //arr->GetTree()->Delete();
521   //delete arr;
522 }
523
524 ///////////////////////////////////////// Point IO  
525 Bool_t AliL3FileHandler::AliPoints2Binary(){
526   Bool_t out = kTRUE;
527   UInt_t npoint;
528   AliL3SpacePointData *data = AliPoints2Memory(npoint);
529   out = Memory2Binary(npoint,data);
530   Free();
531   return out;
532 }
533
534 AliL3SpacePointData * AliL3FileHandler::AliPoints2Memory(UInt_t & npoint){
535   AliL3SpacePointData *data = 0;
536   npoint=0;
537   if(!fInAli){
538     LOG(AliL3Log::kWarning,"AliL3FileHandler::AliPoints2Memory","File")
539     <<"No Input avalible: no object TFile"<<ENDLOG;
540     return 0;
541   }
542   if(!fInAli->IsOpen()){
543     LOG(AliL3Log::kWarning,"AliL3FileHandler::AliPoints2Memory","File")
544     <<"No Input avalible: TFile not opend"<<ENDLOG;
545     return 0;
546   }
547
548   TDirectory *savedir = gDirectory;
549   fInAli->cd();
550   
551   Char_t cname[100];
552   Int_t eventn = 0;
553   sprintf(cname,"TreeC_TPC_%d",eventn);
554   AliTPCClustersArray carray;
555   carray.Setup(fParam);
556   carray.SetClusterType("AliTPCcluster");
557   Bool_t clusterok = carray.ConnectTree(cname);
558   if(!clusterok) return 0;
559
560   AliTPCClustersRow ** clusterrow = 
561                new AliTPCClustersRow*[ (int)carray.GetTree()->GetEntries()];
562   Int_t *rows = new int[ (int)carray.GetTree()->GetEntries()];
563   Int_t *sects = new int[  (int)carray.GetTree()->GetEntries()];
564   Int_t sum=0;
565
566   Int_t lslice,lrow;
567   for(Int_t i=0; i<carray.GetTree()->GetEntries(); i++){
568     AliSegmentID *s = carray.LoadEntry(i);
569     Int_t sector,row;
570     fParam->AdjustSectorRow(s->GetID(),sector,row);
571     rows[i] = row;
572     sects[i] = sector;
573     clusterrow[i] = 0;
574     AliL3Transform::Sector2Slice(lslice,lrow,sector,row);
575     if(fSlice != lslice || lrow<fRowMin || lrow>fRowMax) continue;
576     clusterrow[i] = carray.GetRow(sector,row);
577     if(clusterrow[i])
578       sum+=clusterrow[i]->GetArray()->GetEntriesFast();
579   }
580   UInt_t size = sum*sizeof(AliL3SpacePointData);
581
582   LOG(AliL3Log::kDebug,"AliL3FileHandler::AliPoints2Memory","File")
583   <<AliL3Log::kDec<<"Found "<<sum<<" SpacePoints"<<ENDLOG;
584
585   data = (AliL3SpacePointData *) Allocate(size);
586   npoint = sum;
587   UInt_t n=0; 
588   for(Int_t i=0; i<carray.GetTree()->GetEntries(); i++){
589     if(!clusterrow[i]) continue;
590     Int_t row = rows[i];
591     Int_t sector = sects[i];
592     AliL3Transform::Sector2Slice(lslice,lrow,sector,row);
593     Int_t entries_in_row = clusterrow[i]->GetArray()->GetEntriesFast();
594     for(Int_t j = 0;j<entries_in_row;j++){
595       AliTPCcluster *c = (AliTPCcluster*)(*clusterrow[i])[j];
596       data[n].fZ = c->GetZ();
597       data[n].fY = c->GetY();
598       data[n].fX = fParam->GetPadRowRadii(sector,row);
599       data[n].fID = n+((fSlice&0x7f)<<25)+((fPatch&0x7)<<22);//uli
600       data[n].fPadRow = lrow;
601       data[n].fXYErr = c->GetSigmaY2();
602       data[n].fZErr = c->GetSigmaZ2();
603       if(fMC) fprintf(fMC,"%d %d\n",data[n].fID,c->GetLabel(0));
604       n++;
605     }
606   }
607   for(Int_t i=0;i<carray.GetTree()->GetEntries();i++){
608     Int_t row = rows[i];
609     Int_t sector = sects[i];
610     if(carray.GetRow(sector,row))
611       carray.ClearRow(sector,row);
612   }
613
614   delete [] clusterrow;
615   delete [] rows;
616   delete [] sects;
617   savedir->cd();   
618
619   return data;
620 }
621