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 //-----------------------------------------------------//
20 // Source File : PMDClusteringV2.cxx //
22 // clustering code for alice pmd //
24 //-----------------------------------------------------//
26 /* --------------------------------------------------------------------
27 Code developed by S. C. Phatak, Institute of Physics,
28 Bhubaneswar 751 005 ( phatak@iopb.res.in ) Given the energy deposited
29 ( or ADC value ) in each cell of supermodule ( pmd or cpv ), the code
30 builds up superclusters and breaks them into clusters. The input is
31 in TObjarray and cluster information is in TObjArray.
32 integer clno gives total number of clusters in the supermodule.
33 fClusters is the global ( public ) variables.
34 Others are local ( private ) to the code.
35 At the moment, the data is read for whole detector ( all supermodules
36 and pmd as well as cpv. This will have to be modify later )
37 LAST UPDATE : October 23, 2002
38 -----------------------------------------------------------------------*/
40 #include <Riostream.h>
42 #include <TObjArray.h>
45 #include "AliPMDcludata.h"
46 #include "AliPMDcluster.h"
47 #include "AliPMDisocell.h"
48 #include "AliPMDClustering.h"
49 #include "AliPMDClusteringV2.h"
52 ClassImp(AliPMDClusteringV2)
54 const Double_t AliPMDClusteringV2::fgkSqroot3by2=0.8660254; // sqrt(3.)/2.
56 AliPMDClusteringV2::AliPMDClusteringV2():
57 fPMDclucont(new TObjArray()),
60 for(int i = 0; i < kNDIMX; i++)
62 for(int j = 0; j < kNDIMY; j++)
64 fCoord[0][i][j] = i+j/2.;
65 fCoord[1][i][j] = fgkSqroot3by2*j;
69 // ------------------------------------------------------------------------ //
72 AliPMDClusteringV2::AliPMDClusteringV2(const AliPMDClusteringV2& pmdclv2):
73 AliPMDClustering(pmdclv2),
78 AliError("Copy constructor not allowed ");
81 // ------------------------------------------------------------------------ //
82 AliPMDClusteringV2 &AliPMDClusteringV2::operator=(const AliPMDClusteringV2& /*pmdclv2*/)
85 AliError("Assignment operator not allowed ");
88 // ------------------------------------------------------------------------ //
89 AliPMDClusteringV2::~AliPMDClusteringV2()
93 // ------------------------------------------------------------------------ //
95 void AliPMDClusteringV2::DoClust(Int_t idet, Int_t ismn,
96 Int_t celltrack[48][96],
97 Int_t cellpid[48][96],
98 Double_t celladc[48][96],
99 TObjArray *pmdisocell, TObjArray *pmdcont)
101 // main function to call other necessary functions to do clustering
103 AliPMDcluster *pmdcl = 0;
105 const Float_t ktwobysqrt3 = 1.1547; // 2./sqrt(3.)
106 const Int_t kNmaxCell = 19; // # of cells surrounding a cluster center
107 Int_t i, j, nmx1, incr, id, jd;
110 Int_t celldataX[kNmaxCell], celldataY[kNmaxCell];
111 Int_t celldataTr[kNmaxCell], celldataPid[kNmaxCell];
112 Float_t celldataAdc[kNmaxCell];
114 Double_t cutoff, ave;
115 Double_t edepcell[kNMX];
118 // call the isolated cell search method
120 CalculateIsoCell(idet, ismn, celladc, pmdisocell);
129 else if (ismn >= 12 && ismn <= 23)
135 for (i =0; i < kNMX; i++)
140 for (id = 0; id < ndimXr; id++)
142 for (jd = 0; jd < ndimYr; jd++)
145 i = id + (ndimYr/2-1) - (jd/2);
146 Int_t ij = i + j*kNDIMX;
149 edepcell[ij] = celladc[jd][id];
151 else if (ismn >= 12 && ismn <= 23)
153 edepcell[ij] = celladc[id][jd];
160 TMath::Sort((Int_t)kNMX,edepcell,iord1);// order the data
161 cutoff = fCutoff; // cutoff used to discard cells having ener. dep.
165 for(i = 0;i < kNMX; i++)
171 if(edepcell[i] > cutoff )
177 AliDebug(1,Form("Number of cells having energy >= %f are %d",cutoff,nmx1));
185 AliDebug(1,Form("Number of cells in a SuperM = %d and Average = %f",
188 incr = CrClust(ave, cutoff, nmx1,iord1, edepcell);
189 RefClust(incr,edepcell );
191 Int_t nentries1 = fPMDclucont->GetEntries();
192 AliDebug(1,Form("Detector Plane = %d Serial Module No = %d Number of clusters = %d",idet, ismn, nentries1));
193 AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
194 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
196 AliPMDcludata *pmdcludata =
197 (AliPMDcludata*)fPMDclucont->UncheckedAt(ient1);
198 Float_t cluXC = pmdcludata->GetClusX();
199 Float_t cluYC = pmdcludata->GetClusY();
200 Float_t cluADC = pmdcludata->GetClusADC();
201 Float_t cluCELLS = pmdcludata->GetClusCells();
202 Float_t cluSIGX = pmdcludata->GetClusSigmaX();
203 Float_t cluSIGY = pmdcludata->GetClusSigmaY();
205 Float_t cluY0 = ktwobysqrt3*cluYC;
206 Float_t cluX0 = cluXC - cluY0/2.;
209 // Cluster X centroid is back transformed
213 clusdata[0] = cluX0 - (24-1) + cluY0/2.;
215 else if (ismn >= 12 && ismn <= 23)
217 clusdata[0] = cluX0 - (48-1) + cluY0/2.;
221 clusdata[2] = cluADC;
222 clusdata[3] = cluCELLS;
223 clusdata[4] = cluSIGX;
224 clusdata[5] = cluSIGY;
226 // Cells associated with a cluster
228 for (Int_t ihit = 0; ihit < kNmaxCell; ihit++)
230 Int_t dummyXY = pmdcludata->GetCellXY(ihit);
232 Int_t celldumY = dummyXY%10000;
233 Int_t celldumX = dummyXY/10000;
234 Float_t cellY = (Float_t) celldumY/10;
235 Float_t cellX = (Float_t) celldumX/10;
238 // Cell X centroid is back transformed
242 celldataX[ihit] = (Int_t) ((cellX - (24-1) + cellY/2.) + 0.5);
244 else if (ismn >= 12 && ismn <= 23)
246 celldataX[ihit] = (Int_t) ((cellX - (48-1) + cellY/2.) + 0.5 );
248 celldataY[ihit] = (Int_t) (cellY + 0.5);
250 Int_t irow = celldataX[ihit];
251 Int_t icol = celldataY[ihit];
253 if ((irow >= 0 && irow < 48) && (icol >= 0 && icol < 96))
255 celldataTr[ihit] = celltrack[irow][icol];
256 celldataPid[ihit] = cellpid[irow][icol];
257 celldataAdc[ihit] = (Float_t) celladc[irow][icol];
261 celldataTr[ihit] = -1;
262 celldataPid[ihit] = -1;
263 celldataAdc[ihit] = -1;
268 pmdcl = new AliPMDcluster(idet, ismn, clusdata, celldataX, celldataY,
269 celldataTr, celldataPid, celldataAdc);
272 fPMDclucont->Delete();
274 // ------------------------------------------------------------------------ //
275 Int_t AliPMDClusteringV2::CrClust(Double_t ave, Double_t cutoff, Int_t nmx1,
276 Int_t iord1[], Double_t edepcell[])
278 // Does crude clustering
279 // Finds out only the big patch by just searching the
283 Int_t i,j,k,id1,id2,icl, numcell;
284 Int_t jd1,jd2, icell, cellcount;
285 Int_t clust[2][5000];
286 static Int_t neibx[6] = {1,0,-1,-1,0,1}, neiby[6] = {0,1,1,0,-1,-1};
288 // neibx and neiby define ( incremental ) (i,j) for the neighbours of a
289 // cell. There are six neighbours.
290 // cellcount --- total number of cells having nonzero ener dep
291 // numcell --- number of cells in a given supercluster
293 AliDebug(1,Form("kNMX = %d nmx1 = %d kNDIMX = %d kNDIMY = %d ave = %f cutoff = %f",kNMX,nmx1,kNDIMX,kNDIMY,ave,cutoff));
295 for (j=0; j < kNDIMX; j++)
297 for(k=0; k < kNDIMY; k++)
299 fInfocl[0][j][k] = 0;
300 fInfocl[1][j][k] = 0;
304 for(i=0; i < kNMX; i++)
312 if(edepcell[j] <= cutoff)
314 fInfocl[0][id1][id2] = -1;
317 // ---------------------------------------------------------------
318 // crude clustering begins. Start with cell having largest adc
319 // count and loop over the cells in descending order of adc count
320 // ---------------------------------------------------------------
323 for(icell=0; icell <= nmx1; icell++)
328 if(fInfocl[0][id1][id2] == 0 )
330 // ---------------------------------------------------------------
331 // icl -- cluster #, numcell -- # of cells in it, clust -- stores
332 // coordinates of the cells in a cluster, fInfocl[0][i1][i2] is 1 for
333 // primary and 2 for secondary cells,
334 // fInfocl[1][i1][i2] stores cluster #
335 // ---------------------------------------------------------------
339 fInfocl[0][id1][id2] = 1;
340 fInfocl[1][id1][id2] = icl;
341 fInfcl[0][cellcount] = icl;
342 fInfcl[1][cellcount] = id1;
343 fInfcl[2][cellcount] = id2;
345 clust[0][numcell] = id1;
346 clust[1][numcell] = id2;
347 for(i = 1; i < 5000; i++)
351 // ---------------------------------------------------------------
352 // check for adc count in neib. cells. If ne 0 put it in this clust
353 // ---------------------------------------------------------------
354 for(i = 0; i < 6; i++)
356 jd1 = id1 + neibx[i];
357 jd2 = id2 + neiby[i];
358 if( (jd1 >= 0 && jd1 < kNDIMX) && (jd2 >= 0 && jd2 < kNDIMY) &&
359 fInfocl[0][jd1][jd2] == 0)
362 fInfocl[0][jd1][jd2] = 2;
363 fInfocl[1][jd1][jd2] = icl;
364 clust[0][numcell] = jd1;
365 clust[1][numcell] = jd2;
367 fInfcl[0][cellcount] = icl;
368 fInfcl[1][cellcount] = jd1;
369 fInfcl[2][cellcount] = jd2;
372 // ---------------------------------------------------------------
373 // check adc count for neighbour's neighbours recursively and
374 // if nonzero, add these to the cluster.
375 // ---------------------------------------------------------------
376 for(i = 1;i < 5000; i++)
378 if(clust[0][i] != -1)
382 for(j = 0; j < 6 ; j++)
384 jd1 = id1 + neibx[j];
385 jd2 = id2 + neiby[j];
386 if( (jd1 >= 0 && jd1 < kNDIMX) &&
387 (jd2 >= 0 && jd2 < kNDIMY)
388 && fInfocl[0][jd1][jd2] == 0 )
390 fInfocl[0][jd1][jd2] = 2;
391 fInfocl[1][jd1][jd2] = icl;
393 clust[0][numcell] = jd1;
394 clust[1][numcell] = jd2;
396 fInfcl[0][cellcount] = icl;
397 fInfcl[1][cellcount] = jd1;
398 fInfcl[2][cellcount] = jd2;
407 // ------------------------------------------------------------------------ //
408 void AliPMDClusteringV2::RefClust(Int_t incr, Double_t edepcell[])
410 // Does the refining of clusters
411 // Takes the big patch and does gaussian fitting and
412 // finds out the more refined clusters
414 const Float_t ktwobysqrt3 = 1.1547;
415 const Int_t kNmaxCell = 19;
417 AliPMDcludata *pmdcludata = 0;
420 Int_t i, j, k, i1, i2, id, icl, itest, ihld;
421 Int_t ig, nsupcl, clno, clX,clY;
422 Int_t clxy[kNmaxCell];
425 Double_t x1, y1, z1, x2, y2, z2, rr;
427 Int_t kndim = incr + 1;
434 Double_t *x, *y, *z, *xc, *yc, *zc, *cells, *rcl, *rcs;
436 ncl = new Int_t [kndim];
437 iord = new Int_t [kndim];
438 x = new Double_t [kndim];
439 y = new Double_t [kndim];
440 z = new Double_t [kndim];
441 xc = new Double_t [kndim];
442 yc = new Double_t [kndim];
443 zc = new Double_t [kndim];
444 cells = new Double_t [kndim];
445 rcl = new Double_t [kndim];
446 rcs = new Double_t [kndim];
448 for(Int_t kk = 0; kk < 15; kk++)
450 if( kk < 6 )clusdata[kk] = 0.;
453 // nsupcl = # of superclusters; ncl[i]= # of cells in supercluster i
454 // x, y and z store (x,y) coordinates of and energy deposited in a cell
455 // xc, yc store (x,y) coordinates of the cluster center
456 // zc stores the energy deposited in a cluster, rc is cluster radius
461 for(i = 0; i < kndim; i++)
465 for(i = 0; i <= incr; i++)
467 if(fInfcl[0][i] != nsupcl)
473 AliWarning("RefClust: Too many superclusters!");
480 AliDebug(1,Form("Number of cells = %d Number of Superclusters = %d",
485 for(i = 0; i <= nsupcl; i++)
491 // one cell super-clusters --> single cluster
492 // cluster center at the centyer of the cell
493 // cluster radius = half cell dimension
496 AliWarning("RefClust: Too many clusters! more than 5000");
502 i12 = i1 + i2*kNDIMX;
503 clusdata[0] = fCoord[0][i1][i2];
504 clusdata[1] = fCoord[1][i1][i2];
505 clusdata[2] = edepcell[i12];
512 clY = (Int_t)((ktwobysqrt3*fCoord[1][i1][i2])*10);
513 clX = (Int_t)((fCoord[0][i1][i2] - clY/20.)*10);
514 clxy[0] = clX*10000 + clY ;
516 for(Int_t icltr = 1; icltr < kNmaxCell; icltr++)
520 pmdcludata = new AliPMDcludata(clusdata,clxy);
521 fPMDclucont->Add(pmdcludata);
527 // two cell super-cluster --> single cluster
528 // cluster center is at ener. dep.-weighted mean of two cells
529 // cluster radius == half cell dimension
534 AliWarning("RefClust: Too many clusters! more than 5000");
540 i12 = i1 + i2*kNDIMX;
542 x1 = fCoord[0][i1][i2];
543 y1 = fCoord[1][i1][i2];
549 i12 = i1 + i2*kNDIMX;
551 x2 = fCoord[0][i1][i2];
552 y2 = fCoord[1][i1][i2];
555 clusdata[0] = (x1*z1+x2*z2)/(z1+z2);
556 clusdata[1] = (y1*z1+y2*z2)/(z1+z2);
559 clusdata[4] = (TMath::Sqrt(z1*z2))/(z1+z2);
562 clY = (Int_t)((ktwobysqrt3*y1)*10);
563 clX = (Int_t)((x1 - clY/20.)*10);
564 clxy[0] = clX*10000 + clY ;
566 clY = (Int_t)((ktwobysqrt3*y2)*10);
567 clX = (Int_t)((x2 - clY/20.)*10);
568 clxy[1] = clX*10000 + clY ;
570 for(Int_t icltr = 2; icltr < kNmaxCell; icltr++)
574 pmdcludata = new AliPMDcludata(clusdata, clxy);
575 fPMDclucont->Add(pmdcludata);
580 // super-cluster of more than two cells - broken up into smaller
581 // clusters gaussian centers computed. (peaks separated by > 1 cell)
582 // Begin from cell having largest energy deposited This is first
584 // *****************************************************************
585 // NOTE --- POSSIBLE MODIFICATION: ONE MAY NOT BREAKING SUPERCLUSTERS
586 // IF NO. OF CELLS IS NOT TOO LARGE ( SAY 5 OR 6 )
587 // SINCE WE EXPECT THE SUPERCLUSTER
588 // TO BE A SINGLE CLUSTER
589 //*******************************************************************
593 i12 = i1 + i2*kNDIMX;
595 x[0] = fCoord[0][i1][i2];
596 y[0] = fCoord[1][i1][i2];
597 z[0] = edepcell[i12];
600 for(j = 1; j <= ncl[i]; j++)
606 i12 = i1 + i2*kNDIMX;
608 x[j] = fCoord[0][i1][i2];
609 y[j] = fCoord[1][i1][i2];
610 z[j] = edepcell[i12];
613 // arranging cells within supercluster in decreasing order
614 for(j = 1; j <= ncl[i];j++)
618 for(i1 = 0; i1 < j; i1++)
620 if(itest == 0 && z[iord[i1]] < z[ihld])
623 for(i2 = j-1;i2 >= i1;i2--)
625 iord[i2+1] = iord[i2];
633 // compute the number of clusters and their centers ( first
635 // centers must be separated by cells having smaller ener. dep.
636 // neighbouring centers should be either strong or well-separated
641 for(j = 1; j <= ncl[i]; j++)
646 for(k = 0; k <= ig; k++)
650 rr = Distance(x1,y1,x2,y2);
651 //************************************************************
652 // finetuning cluster splitting
653 // the numbers zc/4 and zc/10 may need to be changed.
654 // Also one may need to add one more layer because our
655 // cells are smaller in absolute scale
656 //************************************************************
659 if( rr >= 1.1 && rr < 1.8 && z[iord[j]] > zc[k]/4.) itest++;
660 if( rr >= 1.8 && rr < 2.1 && z[iord[j]] > zc[k]/10.) itest++;
661 if( rr >= 2.1)itest++;
672 ClustDetails(ncl[i], ig, x, y ,z, xc, yc, zc, rcl, rcs, cells,
676 for(j = 0; j <= ig; j++)
681 AliWarning("RefClust: Too many clusters! more than 5000");
687 clusdata[4] = rcl[j];
688 clusdata[5] = rcs[j];
691 clusdata[3] = ncl[i] + 1;
695 clusdata[3] = cells[j];
698 Int_t ncellcls = testncl[j];
699 if( ncellcls < kNmaxCell )
701 for(Int_t kk = 1; kk <= ncellcls; kk++)
703 Int_t ll = testindex[pp];
704 clY = (Int_t)((ktwobysqrt3*y[ll])*10);
705 clX = (Int_t)((x[ll] - clY/20.)*10);
706 clxy[kk-1] = clX*10000 + clY ;
710 for(Int_t icltr = ncellcls ; icltr < kNmaxCell; icltr++)
715 pmdcludata = new AliPMDcludata(clusdata, clxy);
716 fPMDclucont->Add(pmdcludata);
734 // ------------------------------------------------------------------------ //
735 void AliPMDClusteringV2::ClustDetails(Int_t ncell, Int_t nclust, Double_t x[],
736 Double_t y[], Double_t z[],Double_t xc[],
737 Double_t yc[], Double_t zc[],
738 Double_t rcl[], Double_t rcs[],
739 Double_t cells[], TArrayI &testncl,
745 Int_t kndim1 = ncell + 1;//ncell
747 Int_t kndim3 = nclust + 1;//nclust
749 Int_t i, j, k, i1, i2;
750 Double_t x1, y1, x2, y2, rr, b, c, r1, r2;
751 Double_t sumx, sumy, sumxy, sumxx, sum, sum1, sumyy;
753 Double_t *str, *str1, *xcl, *ycl, *cln;
756 Double_t **clustcell;
757 str = new Double_t [kndim3];
758 str1 = new Double_t [kndim3];
759 xcl = new Double_t [kndim3];
760 ycl = new Double_t [kndim3];
761 cln = new Double_t [kndim3];
763 clustcell = new Double_t *[kndim3];
764 cell = new Int_t *[kndim3];
765 cluster = new Int_t *[kndim1];
766 for(i = 0; i < kndim1; i++)
768 cluster[i] = new Int_t [kndim2];
771 for(i = 0; i < kndim3; i++)
779 cell[i] = new Int_t [kndim2];
780 clustcell[i] = new Double_t [kndim1];
781 for(j = 0; j < kndim1; j++)
785 for(j = 0; j < kndim2; j++)
794 // more than one cluster
795 // checking cells shared between several clusters.
796 // First check if the cell is within
797 // one cell unit ( nearest neighbour). Else,
798 // if it is within 1.74 cell units ( next nearest )
799 // Else if it is upto 2 cell units etc.
801 for (i = 0; i <= ncell; i++)
807 // distance <= 1 cell unit
809 for(j = 0; j <= nclust; j++)
813 rr = Distance(x1, y1, x2, y2);
821 // next nearest neighbour
822 if(cluster[i][0] == 0)
824 for(j=0; j<=nclust; j++)
828 rr = Distance(x1, y1, x2, y2);
829 if(rr <= TMath::Sqrt(3.))
837 // next-to-next nearest neighbour
838 if(cluster[i][0] == 0)
840 for(j=0; j<=nclust; j++)
844 rr = Distance(x1, y1, x2, y2);
854 if(cluster[i][0] == 0)
856 for(j = 0; j <= nclust; j++)
860 rr = Distance(x1, y1, x2, y2);
871 // computing cluster strength. Some cells are shared.
872 for(i = 0; i <= ncell; i++)
874 if(cluster[i][0] != 0)
877 for(j = 1; j <= i1; j++)
885 for(k = 0; k < 5; k++)
887 for(i = 0; i <= ncell; i++)
889 if(cluster[i][0] != 0)
893 for(j = 1; j <= i1; j++)
895 sum += str[cluster[i][j]];
898 for(j = 1; j <= i1; j++)
901 str1[i2] += z[i]*str[i2]/sum;
902 clustcell[i2][i] = z[i]*str[i2]/sum;
908 for(j = 0; j <= nclust; j++)
915 for(i = 0; i <= nclust; i++)
921 for(j = 0; j <= ncell; j++)
923 if(clustcell[i][j] != 0)
925 sumx += clustcell[i][j]*x[j];
926 sumy += clustcell[i][j]*y[j];
927 sum += clustcell[i][j];
928 sum1 += clustcell[i][j]/z[j];
931 //** xcl and ycl are cluster centroid positions ( center of gravity )
939 for(j = 0; j <= ncell; j++)
941 sumxx += clustcell[i][j]*(x[j]-xcl[i])*(x[j]-xcl[i])/sum;
942 sumyy += clustcell[i][j]*(y[j]-ycl[i])*(y[j]-ycl[i])/sum;
943 sumxy += clustcell[i][j]*(x[j]-xcl[i])*(y[j]-ycl[i])/sum;
946 c = sumxx*sumyy-sumxy*sumxy;
947 // ******************r1 and r2 are major and minor axes ( r1 > r2 ).
948 r1 = b/2.+TMath::Sqrt(b*b/4.-c);
949 r2 = b/2.-TMath::Sqrt(b*b/4.-c);
950 // final assignments to proper external variables
960 //To get the cell position in a cluster
962 for(Int_t ii=0; ii<= ncell; ii++)
964 Int_t jj = cluster[ii][0];
965 for(Int_t kk=1; kk<= jj; kk++)
967 Int_t ll = cluster[ii][kk];
969 cell[ll][cell[ll][0]] = ii;
973 testncl.Set(nclust+1);
976 for(Int_t ii=0; ii <= nclust; ii++)
978 testncl[ii] = cell[ii][0];
979 counter += testncl[ii];
981 testindex.Set(counter);
983 for(Int_t ii=0; ii<= nclust; ii++)
985 for(Int_t jj = 1; jj<= testncl[ii]; jj++)
987 Int_t kk = cell[ii][jj];
1001 for(j = 0; j <= ncell; j++)
1014 for(j = 0; j <= ncell; j++)
1016 sumxx += clustcell[i][j]*(x[j]-xcl[i])*(x[j]-xcl[i])/sum;
1017 sumyy += clustcell[i][j]*(y[j]-ycl[i])*(y[j]-ycl[i])/sum;
1018 sumxy += clustcell[i][j]*(x[j]-xcl[i])*(y[j]-ycl[i])/sum;
1021 c = sumxx*sumyy-sumxy*sumxy;
1022 r1 = b/2.+ TMath::Sqrt(b*b/4.-c);
1023 r2 = b/2.- TMath::Sqrt(b*b/4.-c);
1025 // To get the cell position in a cluster
1026 testncl.Set(nclust+1);
1027 testindex.Set(ncell+1);
1028 cell[0][0] = ncell + 1;
1029 testncl[0] = cell[0][0];
1031 for(Int_t ii = 1; ii <= ncell; ii++)
1034 Int_t kk = cell[0][ii];
1038 // final assignments
1046 for(i = 0; i < kndim3; i++)
1048 delete [] clustcell[i];
1051 delete [] clustcell;
1053 for(i = 0; i <kndim1 ; i++)
1055 delete [] cluster[i];
1065 // ------------------------------------------------------------------------ //
1066 Double_t AliPMDClusteringV2::Distance(Double_t x1, Double_t y1,
1067 Double_t x2, Double_t y2)
1069 return TMath::Sqrt((x1-x2)*(x1-x2) + (y1-y2)*(y1-y2));
1071 // ------------------------------------------------------------------------ //
1072 void AliPMDClusteringV2::CalculateIsoCell(Int_t idet, Int_t ismn, Double_t celladc[][96], TObjArray *pmdisocell)
1074 // Does isolated cell search for offline calibration
1076 AliPMDisocell *isocell = 0;
1078 const Int_t kMaxRow = 48;
1079 const Int_t kMaxCol = 96;
1080 const Int_t kCellNeighbour = 6;
1084 Int_t neibx[6] = {1,0,-1,-1,0,1};
1085 Int_t neiby[6] = {0,1,1,0,-1,-1};
1088 for(Int_t irow = 0; irow < kMaxRow; irow++)
1090 for(Int_t icol = 0; icol < kMaxCol; icol++)
1092 if(celladc[irow][icol] > 0)
1095 for(Int_t ii = 0; ii < kCellNeighbour; ii++)
1097 id1 = irow + neibx[ii];
1098 jd1 = icol + neiby[ii];
1099 Float_t adc = (Float_t) celladc[id1][jd1];
1103 if(isocount == kCellNeighbour)
1105 Float_t cadc = (Float_t) celladc[irow][icol];
1107 isocell = new AliPMDisocell(idet,ismn,irow,icol,cadc);
1108 pmdisocell->Add(isocell);
1112 } // neigh cell cond.
1119 // ------------------------------------------------------------------------ //
1120 void AliPMDClusteringV2::SetEdepCut(Float_t decut)
1124 // ------------------------------------------------------------------------ //