1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 Revision 1.7 2002/02/05 09:12:26 hristov
19 Small mods for gcc 3.02
21 Revision 1.6 2001/10/21 19:07:24 hristov
22 Several pointers were set to zero in the default constructors to avoid memory management problems
24 Revision 1.5 2001/01/26 19:57:22 hristov
25 Major upgrade of AliRoot code
27 Revision 1.4 2000/10/05 16:08:15 kowal2
28 Changes due to a new class AliComplexCluster. Forward declarations.
30 Revision 1.3 2000/07/10 20:57:39 hristov
31 Update of TPC code and macros by M.Kowalski
33 Revision 1.2 2000/06/30 12:07:49 kowal2
34 Updated from the TPC-PreRelease branch
36 Revision 1.1.2.1 2000/06/25 08:52:51 kowal2
37 replacing AliClusterFinder
41 //-----------------------------------------------------------------------------
43 // Implementation of class ALITPCCLUSTERFINDER
45 //Class for cluster finding in two dimension.
46 //In the present there are implemented two algorithm
47 //primitive recursion algorithm. (FindPeaks)
48 //Algorithm is not working in case of overlaping clusters
49 //Maximum - minimum in direction algoritm (Find clusters)
50 //In this algoritm we suppose that each cluster has local
51 //maximum. From this local maximum I mus see each point
53 //From maximum i can accept every point in radial
54 //direction which is before border in direction
55 //Border in direction occur if we have next in
56 //direction nder threshold or response begin
57 //to increase in given radial direction
58 //-----------------------------------------------------------------------------
61 #include "AliArrayI.h"
62 #include "TClonesArray.h"
67 #include "AliComplexCluster.h"
68 #include "AliTPCClusterFinder.h"
69 #include <Riostream.h>
70 #include <Riostream.h>
72 //direction constants possible direction in 8 different sectors
76 const Int_t kClStackSize =1000;
81 static AliTPCClusterFinder * gClusterFinder; //for fitting routine
83 void gauss(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
85 AliArrayI * points = gClusterFinder->GetStack();
86 const Int_t nbins = gClusterFinder->GetStackIndex();
91 for (i=0;i<nbins; i++) {
92 Float_t x = points->At(i*3);
93 Float_t y = points->At(i*3+1);
94 Float_t z = points->At(i*3+2);
95 Float_t deltax2 = (x-par[1]);
98 Float_t deltay2 = (y-par[2]);
102 delta = z-par[0]*TMath::Exp(-deltax2-deltay2);
103 chisq += delta*delta;
109 ClassImp(AliTPCClusterFinder)
112 AliTPCClusterFinder::AliTPCClusterFinder()
119 fMulSigma2 = 16; //4 sigma
124 fStack = new AliArrayI;
125 fStack->Set(kClStackSize);
135 fMinuit= new TMinuit(5);
136 fMinuit->SetFCN(gauss);
137 gClusterFinder = this;
142 AliTPCClusterFinder::~AliTPCClusterFinder()
144 if (fDigits != 0) delete fDigits;
147 void AliTPCClusterFinder::SetSigmaX(Float_t s0, Float_t s1x, Float_t s1y)
154 void AliTPCClusterFinder::SetSigmaY(Float_t s0, Float_t s1x, Float_t s1y)
163 Bool_t AliTPCClusterFinder::SetSigma2(Int_t i, Int_t j, Float_t & sigmax2, Float_t &sigmay2)
166 //set sigmax2 and sigma y2 accordig i and j position of cell
169 // Float_t x[3] = {ItoX(i),JtoY(j),0};
173 sigmax2= fSigmaX[0]+fSigmaX[1]*x+fSigmaX[2]*y;
174 sigmay2= fSigmaY[0]+fSigmaY[1]*x+fSigmaY[2]*y;
179 Bool_t AliTPCClusterFinder::SetSigma2(Int_t i, Int_t j, Float_t & sigmax2, Float_t &sigmay2)
182 //set sigmax2 and sigma y2 accordig i and j position of cell
184 if (fDetectorParam==0) {
189 Float_t x[3] = {ItoX(i),JtoY(j),0};
191 fDetectorParam->GetClusterSize(x,fDetectorIndex,0,0,sigma);
192 sigmax2=sigma[0]*(fX2-fX1)*(fX2-fX1)/(fDimX*fDimX);
193 sigmay2=sigma[1]*(fY2-fY1)*(fY2-fY1)/(fDimY*fDimY);
199 void AliTPCClusterFinder::GetHisto(TH2F * his2)
202 UInt_t idim =his2->GetNbinsX();
203 UInt_t jdim =his2->GetNbinsY();
204 fX1 = his2->GetXaxis()->GetXmin();
205 fX2 = his2->GetXaxis()->GetXmax();
206 fY1 = his2->GetYaxis()->GetXmin();
207 fY2 = his2->GetYaxis()->GetXmax();
209 if ( (idim>0) && (jdim>0))
214 Int_t size =idim*jdim;
215 if (fDigits !=0) delete fDigits;
216 fDigits = (Int_t*) new Int_t[size];
217 fCells = (AliCell*) new AliCell[size];
221 for (Int_t i = 0; i<(Int_t)idim;i++)
222 for (Int_t j = 0; j<(Int_t)jdim;j++)
224 Int_t index = his2->GetBin(i+1,j+1);
225 //AliCell * cell = GetCell(i,j);
226 //if (cell!=0) cell->SetSignal(his2->GetBinContent(index));
227 SetSignal(static_cast<int>(his2->GetBinContent(index)),i,j);
235 void AliTPCClusterFinder::FindMaxima()
237 for (Int_t i=0; i<fDimX; i++)
238 for (Int_t j=0;j<fDimY; j++)
239 if (IsMaximum(i,j)) cout<<i<<" "<<j<<"\n";
243 void AliTPCClusterFinder::Transform(AliDigitCluster * c)
245 //transform coordinata from bin coordinata to "normal coordinata"
246 //for example if we initialize finder with histogram
247 //it transform values from bin coordinata to the histogram coordinata
250 c->fMaxX=ItoX(c->fMaxX);
251 c->fMaxY=JtoY(c->fMaxY);
253 c->fSigmaX2=c->fSigmaX2*(fX2-fX1)*(fX2-fX1)/(fDimX*fDimX);
254 c->fSigmaY2=c->fSigmaY2*(fY2-fY1)*(fY2-fY1)/(fDimY*fDimY);
255 c->fArea =c->fArea*(fX2-fX1)*(fY2-fY1)/(fDimX*fDimY);
257 void AliTPCClusterFinder::AddToStack(Int_t i, Int_t j, Int_t signal)
262 if ( ((fStackIndex+2)>=kClStackSize) || (fStackIndex<0) ) return;
263 fStack->AddAt(i,fStackIndex);
264 fStack->AddAt(j,fStackIndex+1);
265 fStack->AddAt(signal,fStackIndex+2);
269 void AliTPCClusterFinder::GetClusterStatistic(AliDigitCluster & cluster)
272 //calculate statistic of cluster
274 Double_t sumxw,sumyw,sumx2w,sumy2w,sumxyw,sumw;
275 Int_t minx,maxx,miny,maxy;
276 sumxw=sumyw=sumx2w=sumy2w=sumxyw=sumw=0;
281 Int_t x0=fStack->At(0);
282 Int_t y0=fStack->At(1);
285 Int_t maxQ=fStack->At(2);
288 for (Int_t i = 0; i<fStackIndex;i+=3){
289 Int_t x = fStack->At(i);
290 Int_t y = fStack->At(i+1);
293 Int_t w = fStack->At(i+2);
312 cluster.fX = sumxw/sumw;
313 cluster.fY = sumyw/sumw;
315 cluster.fSigmaX2 = sumx2w/sumw-cluster.fX*cluster.fX;
316 cluster.fSigmaY2 = sumy2w/sumw-cluster.fY*cluster.fY;
317 cluster.fSigmaXY = sumxyw/sumw-cluster.fX*cluster.fY;
318 cluster.fMaxX = maxQx;
319 cluster.fMaxY = maxQy;
321 cluster.fArea = fStackIndex/3;
322 cluster.fNx = maxx-minx+1;
323 cluster.fNy = maxy-miny+1;
328 void AliTPCClusterFinder::GetClusterFit(AliDigitCluster & cluster)
331 //calculate statistic of cluster
333 Double_t arglist[10];
337 fMinuit->mnexcm("SET ERR", arglist ,1,ierflg);
339 //fistly find starting parameters
340 Int_t minx,maxx,miny,maxy,maxQ,maxQx,maxQy;
342 Float_t sumxw,sumyw,sumw;
352 for (Int_t i = 0; i<fStackIndex;i+=3){
353 Int_t x = fStack->At(i);
354 Int_t y = fStack->At(i+1);
355 Int_t w = fStack->At(i+2);
372 Int_t nx = maxx-minx+1;
373 Int_t ny = maxy-miny+1;
375 SetSigma2(maxQx,maxQy,fCurrentSigmaX2,fCurrentSigmaY2);
376 Double_t vstart[5]={maxQ,sumxw/sumw,sumyw/sumw,1/(2*fCurrentSigmaX2),1/(2*fCurrentSigmaY2)};
377 Double_t step[5]={1.,0.01,0.01,0.01,0.01};
378 fMinuit->mnparm(0, "amp", vstart[0], step[0], 0,0,ierflg);
379 fMinuit->mnparm(1, "x0", vstart[1], step[1], 0,0,ierflg);
380 fMinuit->mnparm(2, "y0", vstart[2], step[2], 0,0,ierflg);
381 fMinuit->mnparm(3, "sx2", vstart[3], step[3], 0,0,ierflg);
382 fMinuit->mnparm(4, "sy2", vstart[4], step[4], 0,0,ierflg);
386 fMinuit->mnfree(0); //set unfixed all parameters
387 //if we have area less then
388 if (over<=21) { //if we dont't have more then 7 points
389 fMinuit->FixParameter(3);
390 fMinuit->FixParameter(4);
393 if (nx<3) fMinuit->FixParameter(3); //fix sigma x if no data in x direction
394 if (ny<3) fMinuit->FixParameter(4); //fix sigma y if no data in y direction
396 fMinuit->mnexcm("MIGRAD", arglist ,2,ierflg);
401 fMinuit->GetParameter(0,x[0],error[0]);
402 fMinuit->GetParameter(1,x[1],error[1]);
403 fMinuit->GetParameter(2,x[2],error[2]);
404 fMinuit->GetParameter(3,x[3],error[3]);
405 fMinuit->GetParameter(4,x[4],error[4]);
409 cluster.fMaxX = maxQx;
410 cluster.fMaxY = maxQy;
413 cluster.fSigmaX2 = 1/TMath::Sqrt(2*x[3]);
414 cluster.fSigmaY2 = 1/TMath::Sqrt(2*x[4]);
415 cluster.fSigmaXY = 0;
417 cluster.fArea = over;
423 Bool_t AliTPCClusterFinder::CheckIfDirBorder(Float_t x, Float_t y,
427 //function which control if given cell with index i, j is the
428 //minimum in direction
429 // x and y are estimate of local maximum
430 //direction is given by the
432 AliCell * cellor= GetCell(i,j);
433 Int_t sigor = GetSignal(i,j);
435 //control derivation in direction
436 //if function grows up in direction then there is border
439 Float_t dd = TMath::Sqrt(dx*dx+dy*dy);
440 Float_t ddx = TMath::Abs(dx);
441 ddx = (ddx>0.5) ? ddx-0.5: 0;
443 Float_t ddy = TMath::Abs(dy);
444 ddy = (ddy>0.5) ? ddy-0.5: 0;
446 Float_t d2 = ddx/(2*fDirSigmaFac*fCurrentSigmaX2)+ddy/(2*fDirSigmaFac*fCurrentSigmaY2); //safety factor
447 //I accept sigmax and sigma y bigge by factor sqrt(fDirsigmaFac)
448 Float_t amp = TMath::Exp(-d2)*fCurrentMaxAmp*fDirAmpFac; //safety factor fDirFac>1
450 if (sigor>amp) return kTRUE;
451 if (dd==0) return kFALSE;
455 virtualcell = GetVirtualSignal(i+dx,j+dy);
456 if (virtualcell <=fThreshold) return kFALSE;
457 if (virtualcell>sigor)
458 if (virtualcell>(sigor+fNoiseTh))
459 {cellor->SetDirBorder(fIndex); return kTRUE;}
462 virtualcell = GetVirtualSignal(i+2*dx,j+2*dy);
463 if (virtualcell>sigor)
464 { cellor->SetDirBorder(fIndex); return kTRUE;}
473 Bool_t AliTPCClusterFinder::IsMaximum(Int_t i, Int_t j)
475 //there is maximum if given digits is 1 sigma over all adjacent
477 //or ther exist virual maximum
478 //is maximum on 24 points neighboring
479 // Bool_t res = kFALSE;
483 AliCell * cell = GetCell(i,j);
484 Int_t signal = GetSignal(i,j);
485 if (cell == 0) return kFALSE;
486 for ( Int_t di=-1;di<=1;di++)
487 for ( Int_t dj=-1;dj<=1;dj++){
488 if ( (di!=0) || (dj!=0))
490 AliCell * cell2=GetCell(i+di,j+dj);
491 Int_t signal2 = GetSignal(i+di,j+dj);
498 if (signal2>signal) return kFALSE;
499 if (signal2>fThreshold) overth++;
500 if (signal2==signal) {
501 if (di<0) return kFALSE;
502 if ( (di+dj)<0) return kFALSE;
504 // if (signal>=signal2){
506 if (signal>fNoiseTh+signal2)
512 //if I have only one neighborough over threshold
513 if (overth<2) return kFALSE;
514 if (over<8) return kFALSE;
518 fCurrentMaxAmp =signal;
519 SetMaximum(fIndex,i,j);
522 //check if there exist virtual maximum
523 for (Float_t ddi=0.;(ddi<1.);ddi+=0.5)
524 for (Float_t ddj=0.;(ddj<1.);ddj+=0.5)
525 if (IsVirtualMaximum(Float_t(i)+ddi,Float_t(j)+ddj)){
526 fCurrentMaxX = i+ddi;
527 fCurrentMaxY = j+ddj;
528 fCurrentMaxAmp =signal;
529 SetMaximum(fIndex,i,j);
535 Bool_t AliTPCClusterFinder::IsVirtualMaximum(Float_t x, Float_t y)
537 //there is maximum if given digits is 1 sigma over all adjacent
539 //is maximum on 24 points neighboring
544 Float_t virtualcell = GetVirtualSignal(x,y);
545 if (virtualcell < 0) return kFALSE;
546 for ( Int_t di=-1;di<=1;di++)
547 for ( Int_t dj=-1;dj<=1;dj++)
548 if ( (di!=0) || (dj!=0))
550 Float_t virtualcell2=GetVirtualSignal(x+di,y+dj);
551 if (virtualcell2 < 0) {
557 if (virtualcell2>fThreshold) overth++;
558 if (virtualcell>=virtualcell2){
560 if (virtualcell>fNoiseTh+virtualcell2)
565 if (overth<2) return kFALSE;
566 //if there exist only one or less neighboring above threshold
567 if (oversigma==8) res = kTRUE;
568 else if ((over==8)&&(GetNType()==8)) res=kTRUE;
570 for ( Int_t di=-2;di<=2;di++)
571 for ( Int_t dj=-2;dj<=2;dj++)
572 if ( (di==2)||(di==-2) || (dj==2)|| (dj==-2) )
574 Float_t virtualcell2=GetVirtualSignal(x+di,y+dj);
575 if (virtualcell2 < 0) {
581 if (virtualcell>=virtualcell2) over+=1;
584 if (over == 24) res=kTRUE;
590 void AliTPCClusterFinder::ResetSignal()
593 Int_t size = fDimX*fDimY;
595 if (rOK==kTRUE) for (Int_t i=0 ; i<size;i++) dig[i] = 0;
600 void AliTPCClusterFinder::ResetStatus()
602 //reset status of signals to not used
603 Int_t size = fDimX*fDimY;
605 if (rOK==kTRUE) for (Int_t i=0 ; i<size;i++)
610 AliCell * AliTPCClusterFinder::GetCell(Int_t i, Int_t j)
612 //return reference to the cell with index i,j
614 if ( (i>=0) && (i<fDimX) && (j>=0) && (j<fDimY) )
615 return &fCells[i+j*fDimX];
619 Float_t AliTPCClusterFinder::GetVirtualSignal(Float_t ri, Float_t rj)
621 //it generate virtual cell as mean value from different cels
622 //after using it must be destructed !!!
625 Int_t ddi = (ri>i)? 1:0;
626 Int_t ddj = (rj>j)? 1:0;
629 for (Int_t di=0;di<=ddi;di++)
630 for (Int_t dj=0;dj<=ddj;dj++)
632 Float_t w = (ri-i-di)*(ri-i-di)+(rj-j-dj)*(rj-j-dj);
633 if (w>0) w=1/TMath::Sqrt(w);
635 AliCell * cel2 =GetCell(i+di,j+dj);
636 Int_t signal2 = GetSignal(i+di,j+dj);
642 if (sumw>0) return (sum/sumw);
649 void AliTPCClusterFinder::SetBlockIndex(Int_t * index)
652 //calculate which indexes we must check for border
654 if (TMath::Abs(index[0])<2) index[2] = 0;
656 index[2] = TMath::Abs(index[0])-1;
657 if (index[0]<0) index[2]*=-1; //first x block
659 if (TMath::Abs(index[1])<2) index[3] = 0;
661 index[3] = TMath::Abs(index[1])-1;
662 if (index[1]<0) index[3]*=-1; //first y block
664 if (TMath::Abs(index[0])<TMath::Abs(index[1])){
669 if (index[0]==index[1]) {
680 //***********************************************************************
681 //***********************************************************************
683 TClonesArray * AliTPCClusterFinder::FindPeaks1(TClonesArray *arr)
685 //find peaks and write it in form of AliTPCcluster to array
687 fClustersArray=new TClonesArray("AliDigitCluster",300);
691 fClustersArray = arr;
692 fIndex = fClustersArray->GetEntriesFast();
697 for (Int_t i=0; i<fDimX; i++)
698 for (Int_t j=0;j<fDimY; j++)
702 AliCell * cell = GetCell(i,j);
703 if (!(cell->IsChecked())) Adjacent(i,j);
704 //if there exists more then 2 digits cluster
705 if (fStackIndex >2 ){
706 if (fBFit==kFALSE) GetClusterStatistic(c);
707 else GetClusterFit(c);
708 //write some important chracteristic area of cluster
711 //write cluster information to array
712 TClonesArray &lclusters = *fClustersArray;
713 new (lclusters[fIndex++]) AliDigitCluster(c);
714 // cout<<"fx="<<c.fX<<" fy"<<c.fY<<"\n";
717 return fClustersArray;
721 TClonesArray * AliTPCClusterFinder::FindPeaks2(TClonesArray *arr)
723 //find peaks and write it in form of AliTPCcluster to array
725 fClustersArray=new TClonesArray("AliDigitCluster",300);
729 fClustersArray = arr;
730 fIndex = fClustersArray->GetEntriesFast();
736 for (Int_t i=0; i<fDimX; i++)
737 for (Int_t j=0;j<fDimY; j++)
740 if (IsMaximum(i,j) == kTRUE){
741 SetSigma2(i,j,fCurrentSigmaX2,fCurrentSigmaY2);
744 //if there exists more then 2 digits cluster
745 if (fStackIndex >2 ){
746 if (fBFit==kFALSE) GetClusterStatistic(c);
747 else GetClusterFit(c);
748 //write some important chracteristic area of cluster
751 //write cluster information to array
752 TClonesArray &lclusters = *fClustersArray;
753 new(lclusters[fIndex++]) AliDigitCluster(c);
754 // cout<<"fx="<<c.fX<<" fy"<<c.fY<<"\n";
758 return fClustersArray;
762 TClonesArray * AliTPCClusterFinder::FindPeaks3(TClonesArray *arr)
764 //find peaks and write it in form of AliTPCcluster to array
766 fClustersArray=new TClonesArray("AliDigitCluster",300);
770 fClustersArray = arr;
771 fIndex = fClustersArray->GetEntriesFast();
779 for (Int_t i=0; i<fDimX; i++)
780 for (Int_t j=0;j<fDimY; j++)
783 if (IsMaximum(i,j) == kTRUE){
784 SetSigma2(i,j,fCurrentSigmaX2,fCurrentSigmaY2);
785 AddToStack(i,j,GetSignal(i,j));
787 //loop over different distance
789 for ( Int_t dd =1;((dd<=dmax) && (naccepted>0));dd++){
791 for (Int_t di = -dd;di<=dd;di++){
792 Int_t ddj = dd-TMath::Abs(di);
793 Int_t sigstart = (ddj>0) ? -1 : 0;
794 for (Int_t sig = sigstart;sig<=1;sig+=2){
796 AliCell *cell= GetCell(i+di,j+dj);
797 Int_t signal = GetSignal(i+di,j+dj);
798 if (cell==0) continue;
803 SetBlockIndex(index); //adjust index to control
804 if ( IsBorder(fIndex,i+index[2],j+index[3]) ||
805 IsBorder(fIndex,i+index[4],j+index[5])) {
806 cell->SetBorder(fIndex);
810 if ( signal<=fThreshold ){
812 cell->SetThBorder(fIndex);
813 if (fBFit==kTRUE) AddToStack(i+di,j+dj,signal);
817 if (CheckIfDirBorder(fCurrentMaxX,fCurrentMaxY,i+di,j+dj) == kTRUE) {
818 if (fBFit==kFALSE) AddToStack(i+di,j+dj,signal/2);
821 AddToStack(i+di,j+dj,signal);
827 } //if there is maximum
828 //if there exists more then 2 digits cluster
829 if (fStackIndex >2 ){
830 if (fBFit==kFALSE) GetClusterStatistic(c);
831 else GetClusterFit(c);
832 //write some important chracteristic area of cluster
835 //write cluster information to array
836 TClonesArray &lclusters = *fClustersArray;
837 new(lclusters[fIndex++]) AliDigitCluster(c);
838 // cout<<"fx="<<c.fX<<" fy"<<c.fY<<"\n";
840 } //lopp over all digits
842 return fClustersArray;
850 void AliTPCClusterFinder::Adjacent(Int_t i,Int_t j)
853 //recursive agorithm program
855 if (fBDistType==kTRUE) {
856 Float_t delta = (i-fCurrentMaxX)*(i-fCurrentMaxX)/fCurrentSigmaX2;
857 delta+=(j-fCurrentMaxY)*(j-fCurrentMaxY)/fCurrentSigmaY2;
858 if (delta > fMulSigma2) {
859 SetDirBorder(fIndex,i,j);
863 AliCell *cell = GetCell(i,j);
864 Int_t signal = GetSignal(i,j);
866 cell->SetChecked(fIndex);
867 if ( (q>fThreshold) || (fBFit==kTRUE)) AddToStack(i,j,q);
872 newcel = GetCell(i-1,j);
873 if (newcel !=0) if (!newcel->IsChecked(fIndex) ) Adjacent(i-1,j);
874 newcel = GetCell(i,j-1);
875 if (newcel !=0) if (!newcel->IsChecked(fIndex) ) Adjacent(i,j-1);
876 newcel = GetCell(i+1,j);
877 if (newcel !=0) if (!newcel->IsChecked(fIndex) ) Adjacent(i+1,j);
878 newcel = GetCell(i,j+1);
879 if (newcel !=0) if (!newcel->IsChecked(fIndex) ) Adjacent(i,j+1);
881 else cell->SetThBorder(fIndex);
886 AliH2F * AliTPCClusterFinder::DrawHisto( const char *option=0,
887 Float_t x1, Float_t x2, Float_t y1, Float_t y2)
890 //draw digits in given array
894 sprintf(ch,"Cluster finder digits ");
895 if ( (fDimX<1)|| (fDimY<1)) {
898 AliH2F * his = new AliH2F(ch,ch,fDimX,fX1,fX2,fDimY,fY1,fY2);
899 //set histogram values
900 for (Int_t i = 0; i<fDimX;i++)
901 for (Int_t j = 0; j<fDimY;j++){
904 his->Fill(x,y,GetSignal(i,j));
907 AliH2F *h2fsub = his->GetSubrange2d(x1,x2,y1,y2);
911 if (his==0) return 0;
912 if (option!=0) his->Draw(option);
918 void AliTPCClusterFinder::DrawCluster(
919 Int_t color, Int_t size, Int_t style)
922 if (fClustersArray==0) return;
923 //draw marker for each of cluster
924 Int_t ncl=fClustersArray->GetEntriesFast();
925 for (Int_t i=0;i<ncl;i++){
926 AliComplexCluster *cl = (AliComplexCluster*)fClustersArray->UncheckedAt(i);
927 TMarker * marker = new TMarker;
928 marker->SetX(cl->fX);
929 marker->SetY(cl->fY);
930 marker->SetMarkerSize(size);
931 marker->SetMarkerStyle(style);
932 marker->SetMarkerColor(color);
939 AliH2F * AliTPCClusterFinder::DrawBorders( const char *option, AliH2F *h, Int_t type ,
940 Float_t x1, Float_t x2, Float_t y1, Float_t y2)
943 //draw digits in given array
947 sprintf(ch,"Cluster finder digits borders");
948 if ( (fDimX<1)|| (fDimY<1)) {
953 else his = new AliH2F(ch,ch,fDimX,fX1,fX2,fDimY,fY1,fY2);
954 //set histogram values
955 for (Int_t i = 0; i<fDimX;i++)
956 for (Int_t j = 0; j<fDimY;j++){
959 if (((type==1)||(type==0))&&IsMaximum(0,i,j)) his->Fill(x,y,16);
960 if (((type==3)||(type==0))&&(IsDirBorder(0,i,j))) his->Fill(x,y,8);
961 if (((type==4)||(type==0))&&(IsThBorder(0,i,j))) his->Fill(x,y,4);
962 if (((type==2)||(type==0))&&IsBorder(0,i,j)) his->Fill(x,y,1);
967 AliH2F *h2fsub = his->GetSubrange2d(x1,x2,y1,y2);
971 if (his==0) return 0;
972 if (option!=0) his->Draw(option);