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