Removed obsolete and no longer maintained cluster finders
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinder.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 // Base Class used to find                                                //
18 // the reconstructed points for ITS                                       //
19 // See also AliITSClusterFinderSPD, AliITSClusterFinderSDD,               //
20 // AliITSClusterFinderSDD  AliITSClusterFinderV2                          //
21 ////////////////////////////////////////////////////////////////////////////
22
23 #include "AliRun.h"
24 #include "AliITSClusterFinder.h"
25 #include "AliITSRecPoint.h"
26 #include "AliITSdigit.h"
27 #include "AliITSDetTypeRec.h"
28 #include "AliITSMap.h"
29 #include "AliITSgeomTGeo.h"
30 #include <TParticle.h>
31 #include "AliMC.h"
32
33 ClassImp(AliITSClusterFinder)
34
35 extern AliRun *gAlice;
36
37 //----------------------------------------------------------------------
38 AliITSClusterFinder::AliITSClusterFinder():
39 TObject(),
40 fDebug(0),
41 fModule(0),
42 fDigits(0),
43 fNdigits(0),
44 fDetTypeRec(0),
45 fClusters(0),
46 fNRawClusters(0),
47 fMap(0),
48 fNperMax(0),
49 fDeclusterFlag(0),
50 fClusterSize(0),
51 fNPeaks(-1),
52 fNModules(AliITSgeomTGeo::GetNModules()),
53 fEvent(0){
54     // default cluster finder
55     // Input:
56     //   none.
57     // Output:
58     //   none.
59     // Return:
60     //   A default constructed AliITSCulsterFinder
61 }
62 //----------------------------------------------------------------------
63 AliITSClusterFinder::AliITSClusterFinder(AliITSDetTypeRec* dettyp):
64 TObject(),
65 fDebug(0),
66 fModule(0),
67 fDigits(0),
68 fNdigits(0),
69 fDetTypeRec(dettyp),
70 fClusters(0),
71 fNRawClusters(0),
72 fMap(0),
73 fNperMax(0),
74 fDeclusterFlag(0),
75 fClusterSize(0),
76 fNPeaks(-1),
77 fNModules(AliITSgeomTGeo::GetNModules()),
78 fEvent(0){
79     // default cluster finder
80     // Standard constructor for cluster finder
81     // Input:
82     //   AliITSsegmentation *seg  The segmentation class to be used
83     //   AliITSresponse     *res  The response class to be used
84     // Output:
85     //   none.
86     // Return:
87     //   A Standard constructed AliITSCulsterFinder
88
89     SetNperMax();
90     SetClusterSize();
91     SetDeclusterFlag();
92 }
93 //----------------------------------------------------------------------
94 AliITSClusterFinder::AliITSClusterFinder(AliITSDetTypeRec* dettyp,
95                                          TClonesArray *digits):
96 TObject(),
97 fDebug(0),
98 fModule(0),
99 fDigits(digits),
100 fNdigits(0),
101 fDetTypeRec(dettyp),
102 fClusters(0),
103 fNRawClusters(0),
104 fMap(0),
105 fNperMax(0),
106 fDeclusterFlag(0),
107 fClusterSize(0),
108 fNPeaks(-1),
109 fNModules(AliITSgeomTGeo::GetNModules()),
110 fEvent(0){
111     // default cluster finder
112     // Standard + cluster finder constructor
113     // Input:
114     //   AliITSsegmentation *seg  The segmentation class to be used
115     //   AliITSresponse     *res  The response class to be used
116     //   TClonesArray    *digits  Array of digits to be used
117     // Output:
118     //   none.
119     // Return:
120     //   A Standard constructed AliITSCulsterFinder
121
122     fNdigits = fDigits->GetEntriesFast();
123     SetNperMax();
124     SetClusterSize();
125     SetDeclusterFlag();
126 }
127
128 //______________________________________________________________________
129 AliITSClusterFinder::AliITSClusterFinder(const AliITSClusterFinder &source) : TObject(source),
130 fDebug(source.fDebug),
131 fModule(source.fModule),
132 fDigits(),
133 fNdigits(source.fNdigits),
134 fDetTypeRec(),
135 fClusters(),
136 fNRawClusters(source.fNRawClusters),
137 fMap(),
138 fNperMax(source.fNperMax),
139 fDeclusterFlag(source.fDeclusterFlag),
140 fClusterSize(source.fClusterSize),
141 fNPeaks(source.fNPeaks),
142 fNModules(source.fNModules),
143 fEvent(source.fEvent) {
144   // Copy constructor
145   // Copies are not allowed. The method is protected to avoid misuse.
146   AliError("Copy constructor not allowed\n");
147 }
148
149
150 //______________________________________________________________________
151 //AliITSClusterFinder& AliITSClusterFinder::operator=(const AliITSClusterFinder& /* source */){
152   // Assignment operator
153   // Assignment is not allowed. The method is protected to avoid misuse.
154 //  Fatal("= operator","Assignment operator not allowed\n");
155 //  return *this;
156 //}
157
158 //----------------------------------------------------------------------
159 AliITSClusterFinder::~AliITSClusterFinder(){
160     // destructor cluster finder
161     // Input:
162     //   none.
163     // Output:
164     //   none.
165     // Return:
166     //   none.
167
168     if(fMap) {delete fMap;}
169     // Zero local pointers. Other classes own these pointers.
170     fMap          = 0;
171     fDigits       = 0;
172     fNdigits      = 0;
173     fNRawClusters = 0;
174     fNperMax      = 0;
175     fDeclusterFlag= 0;
176     fClusterSize  = 0;
177     fNPeaks       = 0;
178     fDetTypeRec   = 0;
179
180 }
181 //__________________________________________________________________________
182 void AliITSClusterFinder::InitGeometry(){
183  //
184  // Initialisation of ITS geometry
185  //
186   Int_t mmax=AliITSgeomTGeo::GetNModules();
187   for (Int_t m=0; m<mmax; m++) {
188     Int_t lay,lad,det; AliITSgeomTGeo::GetModuleId(m,lay,lad,det);
189     fNdet[m] = (lad-1)*AliITSgeomTGeo::GetNDetectors(lay) + (det-1);
190     fNlayer[m] = lay-1;
191   }
192 }
193
194
195
196 //----------------------------------------------------------------------
197 void AliITSClusterFinder::AddCluster(Int_t branch, AliITSRawCluster *c){
198     // Add a raw cluster copy to the list
199     // Input:
200     //   Int_t       branch  The branch to which the cluster is to be added to
201     //   AliITSRawCluster *c The cluster to be added to the array of clusters
202     // Output:
203     //   none.
204     // Return:
205     //   none.
206
207    if(!fDetTypeRec) {
208     Error("AddCluster","fDetTypeRec is null!");
209     return;
210   }
211   fDetTypeRec->AddCluster(branch,c); 
212   fNRawClusters++;
213 }
214 //----------------------------------------------------------------------
215 void AliITSClusterFinder::AddCluster(Int_t branch, AliITSRawCluster *c, 
216                                      AliITSRecPoint &rp){
217     // Add a raw cluster copy to the list and the RecPoint
218     // Input:
219     //   Int_t       branch  The branch to which the cluster is to be added to
220     //   AliITSRawCluster *c The cluster to be added to the array of clusters
221     //   AliITSRecPoint  &rp The RecPoint to be added to the array of RecPoints
222     // Output:
223     //   none.
224     // Return:
225     //   none.
226   if(!fDetTypeRec) {
227     Error("AddCluster","fDetTypeRec is null!");
228     return;
229   }
230
231   fDetTypeRec->AddCluster(branch,c); 
232   fNRawClusters++;
233   fDetTypeRec->AddRecPoint(rp); 
234
235 }
236 /*
237 //______________________________________________________________________
238 void AliITSClusterFinder::CheckLabels(Int_t lab[3]) {
239   //------------------------------------------------------------
240   // Tries to find mother's labels
241   //------------------------------------------------------------
242
243   if(lab[0]<0 && lab[1]<0 && lab[2]<0) return; // In case of no labels just exit
244   // Check if simulation
245   AliMC* mc = gAlice->GetMCApp();
246   if(!mc)return;
247
248   Int_t ntracks = mc->GetNtrack();
249   for (Int_t i=0;i<3;i++){
250     Int_t label = lab[i];
251     if (label>=0 && label<ntracks) {
252       TParticle *part=(TParticle*)mc->Particle(label);
253       if (part->P() < 0.005) {
254         Int_t m=part->GetFirstMother();
255         if (m<0) {      
256           continue;
257         }
258         if (part->GetStatusCode()>0) {
259           continue;
260         }
261         lab[i]=m;       
262       }
263     }    
264   }
265   
266 }
267 */
268 //______________________________________________________________________
269 void AliITSClusterFinder::FindRawClusters(Int_t module){
270     // Default Cluster finder.
271     // Input:
272     //   Int_t module   Module number for which culster are to be found.
273     // Output:
274     //   none.
275     // Return:
276     //   none.
277     const Int_t kelms = 10;
278     Int_t ndigits = fDigits->GetEntriesFast();
279     TObjArray *digs = new TObjArray(ndigits);
280     TObjArray *clusts = new TObjArray(ndigits); // max # cluster
281     TObjArray *clust0=0; // A spacific cluster of digits
282     TObjArray *clust1=0; // A spacific cluster of digits
283     AliITSdigit *dig=0; // locat pointer to a digit
284     Int_t i=0,nc=0,j[4],k,k2=0;
285
286     // Copy all digits for this module into a local TObjArray.
287     for(i=0;i<ndigits;i++) digs->AddAt(new AliITSdigit(*(GetDigit(i))),i);
288     digs->Sort();
289     // First digit is a cluster.
290     i  = 0;
291     nc = 0;
292     clusts->AddAt(new TObjArray(kelms),nc);
293     clust0 = (TObjArray*)(clusts->At(nc));
294     clust0->AddAtFree(digs->At(i)); // move owner ship from digs to clusts
295     nc++;
296     for(i=1;i<ndigits;i++){
297         if(IsNeighbor(digs,i,j)){
298             dig = (AliITSdigit*)(digs->At(j[0]));
299             // Add to existing cluster. Find which cluster this digis 
300             for(k=0;k<nc;k++){
301                 clust0 = ((TObjArray*)(clusts->At(k)));
302                 if(clust0->IndexOf(dig)>=0) break;
303             } // end for k
304             if(k>=nc){
305                 Fatal("FindRawClusters","Digit not found as expected");
306             } // end if
307             if(j[1]>=0){
308                 dig = (AliITSdigit*)(digs->At(j[1]));
309                 // Add to existing cluster. Find which cluster this digis 
310                 for(k2=0;k2<nc;k2++){
311                     clust1 = ((TObjArray*)(clusts->At(k2)));
312                     if(clust1->IndexOf(dig)>=0) break;
313                 } // end for k2
314                 if(k2>=nc){
315                     Fatal("FindRawClusters","Digit not found as expected");
316                 } // end if
317             } // end if j[1]>=0
318             // Found cluster with neighboring digits add this one to it.
319             if(clust0==clust1){ // same cluster
320                 clust0->AddAtFree(digs->At(i));
321                 clust0 = 0; // finished with cluster. zero for safty
322                 clust1 = 0; // finished wit hcluster. zero for safty
323             }else{ // two different clusters which need to be merged.
324                 clust0->AddAtFree(digs->At(i)); // Add digit to this cluster.
325                 for(k=0;k<clust1->GetEntriesFast();k++){
326                     // move clust1 into clust0
327                     //move digit to this cluster
328                     clust0->AddAtFree(clust1->At(k));
329                     clust1->AddAt(0,k); // zero this one
330                 } // end for k
331                 delete clust1;
332                 clusts->AddAt(0,k2); // zero array of clusters element clust1
333                 clust0 = 0; // finished with cluster. zero for safty
334                 clust1 = 0; // finished wit hcluster. zero for safty
335             } // end if clust0==clust1
336         }else{// New cluster
337             clusts->AddAt(new TObjArray(kelms),nc);
338             clust0 = ((TObjArray*)(clusts->At(nc)));
339             // move owner ship from digs to clusts
340             clust0->AddAtFree(digs->At(i));
341             clust0 = 0; // finished with cluster. zero for safty
342             nc++;
343         } // End if IsNeighbor
344     } // end for i
345     // There are now nc clusters in clusts. Each element of clust is an
346     // array of digits which are clustered together.
347
348     // For each cluster call detector specific CreateRecPoints
349     for(i=0;i<nc;i++) CreateRecPoints((TObjArray*)(clusts->At(i)),module);
350
351     // clean up at the end.
352     for(i=0;i<nc;i++){ 
353         clust0 =(TObjArray*)(clusts->At(i));
354         // Digits deleted below, so zero this TObjArray
355         for(k=0;k<clust0->GetEntriesFast();k++) clust0->AddAt(0,k);
356         delete clust0; // Delete this TObjArray
357         clusts->AddAt(0,i); // Contents deleted above, so zero it.
358     } // end for i
359     delete clusts; // Delete this TObjArray/
360     // Delete the digits then the TObjArray which containted them.
361     for(i=0;i<ndigits;i++) delete ((AliITSdigit*)(digs->At(i)));
362     delete digs;
363 }
364 //______________________________________________________________________
365 Bool_t AliITSClusterFinder::IsNeighbor(TObjArray *digs,Int_t i,Int_t n[])const{
366     // Locagical function which checks to see if digit i has a neighbor.
367     // If so, then it returns kTRUE and its neighbor index j.
368     // This routine checks if the digits are side by side or one before the
369     // other. Requires that the array of digits be in proper order.
370     // Returns kTRUE in the following cases.
371     //                 ji   0j   if kdiagonal  j0    0i
372     //                 00   0i   if kdiagonal  0i    j0
373     // Inputs:
374     //    TObjArray *digs   Array to search for neighbors in
375     //    Int_t      i      Index of digit for which we are searching for
376     //                      a neighbor of.
377     // Output:
378     //    Int_t      j[4]   Index of one or more of the digits which is a
379     //                      neighbor of digit a index i.
380     // Return:
381     //    Bool_t            kTRUE if a neighbor was found kFALSE otherwise.
382     Int_t ix,jx,iz,jz,j;
383     const Bool_t kdiagonal=kFALSE;
384     Bool_t nei[4];
385
386     // No neighbors found if array empty.
387     if(digs->GetEntriesFast()<=0) return kFALSE;
388     // can not be a digit with first element or elements out or range
389     if(i<=0 || i>=digs->GetEntriesFast()) return kFALSE;
390
391     for(j=0;j<4;j++){n[j] = -1;nei[j]=kFALSE;}
392     ix = ((AliITSdigit*)(digs->At(i)))->GetCoord1();
393     iz = ((AliITSdigit*)(digs->At(i)))->GetCoord2();
394     for(j=0;j<i;j++){
395         jx = ((AliITSdigit*)(digs->At(j)))->GetCoord1();
396         jz = ((AliITSdigit*)(digs->At(j)))->GetCoord2();
397         if(jx+1==ix && jz  ==iz){n[0] = j;nei[0] = kTRUE;}
398         if(jx  ==ix && jz+1==iz){n[1] = j;nei[1] = kTRUE;}
399         if(jx+1==ix && jz+1==iz){n[2] = j;nei[2] = kTRUE;}
400         if(jx+1==ix && jz-1==iz){n[3] = j;nei[3] = kTRUE;}
401     } // end for k
402     if(nei[0]||nei[1]) return kTRUE;
403     if(kdiagonal&&(nei[2]||nei[3])) return kTRUE;
404     // no Neighbors found.
405     return kFALSE;
406 }
407
408 //______________________________________________________________________
409 void AliITSClusterFinder::Print(ostream *os) const{
410     //Standard output format for this class
411     // Inputs:
412     //    ostream *os   Output stream
413     // Output:
414     //    ostream *os   Output stream
415     // Return:
416     //    none.
417
418     *os << fDebug<<",";
419     *os << fModule<<",";
420     *os << fNdigits<<",";
421     *os << fNRawClusters<<",";
422     *os << fNperMax<<",";
423     *os << fDeclusterFlag<<",";
424     *os << fClusterSize<<",";
425     *os << fNPeaks<<endl;
426 }
427 //______________________________________________________________________
428 void AliITSClusterFinder::Read(istream *is)  {
429     //Standard input for this class
430     // Inputs:
431     //    istream *is   Input stream
432     // Output:
433     //    istream *is   Input stream
434     // Return:
435     //    none.
436
437     *is >> fDebug;
438     *is >> fModule;
439     *is >> fNdigits;
440     *is >> fNRawClusters;
441     *is >> fNperMax;
442     *is >> fDeclusterFlag;
443     *is >> fClusterSize;
444     *is >> fNPeaks;
445 }
446 //______________________________________________________________________
447 ostream &operator<<(ostream &os,AliITSClusterFinder &source){
448     // Standard output streaming function.
449     // Inputs:
450     //    ostream             *os     Output stream
451     //    AliITSClusterFinder &source Class to be printed
452     // Output:
453     //    ostream             *os     Output stream
454     // Return:
455     //    none.
456
457     source.Print(&os);
458     return os;
459 }
460 //______________________________________________________________________
461 istream &operator>>(istream &is,AliITSClusterFinder &source){
462     // Standard output streaming function.
463     // Inputs:
464     //    istream              *is      Input stream
465     //     AliITSClusterFinder &source  Class to be read in.
466     // Output:
467     //    istream              *is      Input stream
468     // Return:
469     //    none.
470
471     source.Read(&is);
472     return is;
473 }
474 //______________________________________________________________________
475 void AliITSClusterFinder::CheckLabels2(Int_t lab[10]) {
476   //------------------------------------------------------------
477   // Tries to find mother's labels
478   //------------------------------------------------------------
479   AliRunLoader *rl = AliRunLoader::GetRunLoader();
480   TTree *trK=(TTree*)rl->TreeK();
481
482   if(trK){
483     Int_t nlabels =0; 
484     for (Int_t i=0;i<10;i++) if (lab[i]>=0) nlabels++;
485     if(nlabels == 0) return; // In case of no labels just exit
486
487
488     Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
489
490     for (Int_t i=0;i<10;i++){
491       Int_t label = lab[i];
492       if (label>=0 && label<ntracks) {
493         TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
494
495         if (part->P() < 0.02) {
496           Int_t m=part->GetFirstMother();
497           if (m<0) {    
498             continue;
499           }
500           if (part->GetStatusCode()>0) {
501             continue;
502           }
503           lab[i]=m;       
504         }
505         else
506           if (part->P() < 0.12 && nlabels>3) {
507             lab[i]=-2;
508             nlabels--;
509           } 
510       }
511       else{
512         if ( (label>ntracks||label <0) && nlabels>3) {
513           lab[i]=-2;
514           nlabels--;
515         } 
516       }
517     }  
518     if (nlabels>3){
519       for (Int_t i=0;i<10;i++){
520         if (nlabels>3){
521           Int_t label = lab[i];
522           if (label>=0 && label<ntracks) {
523             TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
524             if (part->P() < 0.1) {
525               lab[i]=-2;
526               nlabels--;
527             }
528           }
529         }
530       }
531     }
532
533     //compress labels -- if multi-times the same
534     Int_t lab2[10];
535     for (Int_t i=0;i<10;i++) lab2[i]=-2;
536     for (Int_t i=0;i<10  ;i++){
537       if (lab[i]<0) continue;
538       for (Int_t j=0;j<10 &&lab2[j]!=lab[i];j++){
539         if (lab2[j]<0) {
540           lab2[j]= lab[i];
541           break;
542         }
543       }
544     }
545     for (Int_t j=0;j<10;j++) lab[j]=lab2[j];
546   
547   }
548 }
549
550 //______________________________________________________________________
551 void AliITSClusterFinder::AddLabel(Int_t lab[10], Int_t label) {
552   //add label to the cluster
553   AliRunLoader *rl = AliRunLoader::GetRunLoader();
554   TTree *trK=(TTree*)rl->TreeK();
555   if(trK){
556     if(label<0) return; // In case of no label just exit
557
558     Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
559     if (label>ntracks) return;
560     for (Int_t i=0;i<10;i++){
561       //    if (label<0) break;
562       if (lab[i]==label) break;
563       if (lab[i]<0) {
564         lab[i]= label;
565         break;
566       }
567     }
568   }
569 }
570
571
572 //______________________________________________________________________
573 void AliITSClusterFinder:: 
574 FindCluster(Int_t k,Int_t maxz,AliBin *bins,Int_t &n,Int_t *idx) {
575   //------------------------------------------------------------
576   // returns an array of indices of digits belonging to the cluster
577   // (needed when the segmentation is not regular) 
578   //------------------------------------------------------------
579   if (n<200) idx[n++]=bins[k].GetIndex();
580   bins[k].Use();
581
582   if (bins[k-maxz].IsNotUsed()) FindCluster(k-maxz,maxz,bins,n,idx);
583   if (bins[k-1   ].IsNotUsed()) FindCluster(k-1   ,maxz,bins,n,idx);
584   if (bins[k+maxz].IsNotUsed()) FindCluster(k+maxz,maxz,bins,n,idx);
585   if (bins[k+1   ].IsNotUsed()) FindCluster(k+1   ,maxz,bins,n,idx);
586   /*
587   if (bins[k-maxz-1].IsNotUsed()) FindCluster(k-maxz-1,maxz,bins,n,idx);
588   if (bins[k-maxz+1].IsNotUsed()) FindCluster(k-maxz+1,maxz,bins,n,idx);
589   if (bins[k+maxz-1].IsNotUsed()) FindCluster(k+maxz-1,maxz,bins,n,idx);
590   if (bins[k+maxz+1].IsNotUsed()) FindCluster(k+maxz+1,maxz,bins,n,idx);
591   */
592 }
593
594 //______________________________________________________________________
595 Bool_t AliITSClusterFinder::IsMaximum(Int_t k,Int_t max,const AliBin *bins) {
596   //------------------------------------------------------------
597   //is this a local maximum ?
598   //------------------------------------------------------------
599   UShort_t q=bins[k].GetQ();
600   if (q==1023) return kFALSE;
601   if (bins[k-max].GetQ() > q) return kFALSE;
602   if (bins[k-1  ].GetQ() > q) return kFALSE; 
603   if (bins[k+max].GetQ() > q) return kFALSE; 
604   if (bins[k+1  ].GetQ() > q) return kFALSE; 
605   if (bins[k-max-1].GetQ() > q) return kFALSE;
606   if (bins[k+max-1].GetQ() > q) return kFALSE; 
607   if (bins[k+max+1].GetQ() > q) return kFALSE; 
608   if (bins[k-max+1].GetQ() > q) return kFALSE;
609   return kTRUE; 
610 }
611
612 //______________________________________________________________________
613 void AliITSClusterFinder::
614 FindPeaks(Int_t k,Int_t max,AliBin *b,Int_t *idx,UInt_t *msk,Int_t& n) {
615   //------------------------------------------------------------
616   //find local maxima
617   //------------------------------------------------------------
618   if (n<31)
619   if (IsMaximum(k,max,b)) {
620     idx[n]=k; msk[n]=(2<<n);
621     n++;
622   }
623   b[k].SetMask(0);
624   if (b[k-max].GetMask()&1) FindPeaks(k-max,max,b,idx,msk,n);
625   if (b[k-1  ].GetMask()&1) FindPeaks(k-1  ,max,b,idx,msk,n);
626   if (b[k+max].GetMask()&1) FindPeaks(k+max,max,b,idx,msk,n);
627   if (b[k+1  ].GetMask()&1) FindPeaks(k+1  ,max,b,idx,msk,n);
628 }
629
630 //______________________________________________________________________
631 void AliITSClusterFinder::
632 MarkPeak(Int_t k, Int_t max, AliBin *bins, UInt_t m) {
633   //------------------------------------------------------------
634   //mark this peak
635   //------------------------------------------------------------
636   UShort_t q=bins[k].GetQ();
637
638   bins[k].SetMask(bins[k].GetMask()|m); 
639
640   if (bins[k-max].GetQ() <= q)
641      if ((bins[k-max].GetMask()&m) == 0) MarkPeak(k-max,max,bins,m);
642   if (bins[k-1  ].GetQ() <= q)
643      if ((bins[k-1  ].GetMask()&m) == 0) MarkPeak(k-1  ,max,bins,m);
644   if (bins[k+max].GetQ() <= q)
645      if ((bins[k+max].GetMask()&m) == 0) MarkPeak(k+max,max,bins,m);
646   if (bins[k+1  ].GetQ() <= q)
647      if ((bins[k+1  ].GetMask()&m) == 0) MarkPeak(k+1  ,max,bins,m);
648 }
649
650 //______________________________________________________________________
651 void AliITSClusterFinder::
652 MakeCluster(Int_t k,Int_t max,AliBin *bins,UInt_t m,AliITSRecPoint &c) {
653   //------------------------------------------------------------
654   //make cluster using digits of this peak
655   //------------------------------------------------------------
656   Float_t q=(Float_t)bins[k].GetQ();
657   Int_t i=k/max, j=k-i*max;
658
659   c.SetQ(c.GetQ()+q);
660   c.SetY(c.GetY()+i*q);
661   c.SetZ(c.GetZ()+j*q);
662   c.SetSigmaY2(c.GetSigmaY2()+i*i*q);
663   c.SetSigmaZ2(c.GetSigmaZ2()+j*j*q);
664
665   bins[k].SetMask(0xFFFFFFFE);
666   
667   if (bins[k-max].GetMask() == m) MakeCluster(k-max,max,bins,m,c);
668   if (bins[k-1  ].GetMask() == m) MakeCluster(k-1  ,max,bins,m,c);
669   if (bins[k+max].GetMask() == m) MakeCluster(k+max,max,bins,m,c);
670   if (bins[k+1  ].GetMask() == m) MakeCluster(k+1  ,max,bins,m,c);
671 }