StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StSvtClusterFinder.cc
1  /***************************************************************************
2  *
3  * $Id: StSvtClusterFinder.cc,v 1.10 2005/07/23 03:37:33 perev Exp $
4  *
5  * Author: Selemon Bekele
6  ***************************************************************************
7  *
8  * Description: SVT Cluster Finder Class
9  *
10  ***************************************************************************
11  *
12  * $Log: StSvtClusterFinder.cc,v $
13  * Revision 1.10 2005/07/23 03:37:33 perev
14  * IdTruth + Cleanup
15  *
16  * Revision 1.9 2003/09/02 17:59:06 perev
17  * gcc 3.2 updates + WarnOff
18  *
19  * Revision 1.8 2003/01/28 20:28:18 munhoz
20  * new filters for clusters
21  *
22  * Revision 1.7 2001/04/25 18:38:28 perev
23  * HPcorrs
24  *
25  * Revision 1.6 2000/11/30 20:43:17 caines
26  * Use database
27  *
28  * Revision 1.5 2000/10/31 16:19:37 caines
29  * Improved speed of zeroing arrays
30  *
31  * Revision 1.4 2000/08/21 13:06:58 bekele
32  * Improved restting arrays, cluster finder
33  * now runs faster.
34  *
35  * $Log: StSvtClusterFinder.cc,v $
36  * Revision 1.10 2005/07/23 03:37:33 perev
37  * IdTruth + Cleanup
38  *
39  * Revision 1.9 2003/09/02 17:59:06 perev
40  * gcc 3.2 updates + WarnOff
41  *
42  * Revision 1.8 2003/01/28 20:28:18 munhoz
43  * new filters for clusters
44  *
45  * Revision 1.7 2001/04/25 18:38:28 perev
46  * HPcorrs
47  *
48  * Revision 1.6 2000/11/30 20:43:17 caines
49  * Use database
50  *
51  * Revision 1.5 2000/10/31 16:19:37 caines
52  * Improved speed of zeroing arrays
53  *
54  * Revision 1.4 2000/08/21 13:06:58 caines
55  * Much improved hit finding and fitting
56  *
57  * Revision 1.3 2000/07/16 22:32:54 caines
58  * Fills spacepoint table correctly
59  *
60  * Revision 1.2 2000/07/13 14:50:49 caines
61  * Improvements on not saving single pixels
62  *
63  * Revision 1.1 2000/07/06 03:50:33 caines
64  * First version of cluster finder and fitter
65  *
66  **************************************************************************/
67 
68 #include <Stiostream.h>
69 #include "StSvtClusterFinder.hh"
70 #include "StSvtClassLibrary/StSvtHybridData.hh"
71 #include "StSequence.hh"
72 
73 //ClassImp(StSvtClusterFinder)
74 
75 StSvtClusterFinder::StSvtClusterFinder()
76 {
77 
78 hybdata = NULL;
79 sequence = NULL;
80 mAnolist = NULL;
81 cluIndex = 0;
82 mNumOfAnodes = 0;
83 mSequence = 0;
84 
85 
86 int i,j;
87 
88 for ( i=0; i<8000; i++)
89  for ( j=0; j<500; j++)
90  mCluster[i][j] = 0;
91 
92 for ( i=0; i<240; i++)
93  {
94  mNumSeq[i] = 0;
95  for ( j=0; j<128; j++)
96  mSeqFlag[i][j] = 0;
97  }
98 
99 for ( i=0; i<500; i++)
100  {mContainer1[i] = 0;
101  mContainer2[i] = 0;}
102 
103  }
104 
105 
106 StSvtClusterFinder::~StSvtClusterFinder()
107 {
108 }
109 
110 void StSvtClusterFinder::setHybridPointer(StSvtHybridData* hdata)
111 {
112  hybdata = hdata;
113  mNumOfAnodes = hdata->getAnodeList(mAnolist);
114 }
115 
116 
117 void StSvtClusterFinder::ClusterFinder()
118  {
119  cluIndex = 0;
120  //if (m_hybIndex==69){
121  //cout << " In clumaker, the number of anodes is " << mNumOfAnodes << endl;
122  //}
123 
124  for(int mAnode = 0; mAnode<mNumOfAnodes; mAnode++)
125  {
126  hybdata->getListSequences(mAnode,mSequence,sequence); //need to decide where to get sequences
127  //cout<<"mSequence = "<<mSequence<<endl;
128 
129  mNumSeq[mAnode] = mSequence;
130  for(int mSeq = 0; mSeq < mSequence; mSeq++)
131  {
132 
133  if(mSeqFlag[mAnode][mSeq] == 0)
134  getClusterMembers(mAnode,mSeq);
135  }
136  }
137 
138  mNumOfClusters = cluIndex;
139 
140  }
141 
142 
143 void StSvtClusterFinder::getClusterMembers(int& mAnode, int &mSeq)
144  {
145  int members, memCount , breakAn = 0, status;
146  int mSeqStart, mSeqLength, mSeqEnd;
147  int newMem = 0;
148  int mRightFlagCounter = 0;
149  int mLeftFlagCounter = 0;
150 
151  ++cluIndex;
152  memCount = 0;
153  mCluster[cluIndex-1][memCount] = 1000*mAnode + mSeq;
154  //cout<<mAnolist[mAnode]<<endl;
155  //cout<<mSeq<<endl;
156  mContainer1[0] = mCluster[cluIndex-1][memCount];
157  mSeqFlag[mAnode][mSeq] = 1;
158  members = 1;
159  memCount = 1;
160 
161 // mContainer1[0] = mCluster[cluIndex-1][memCount];
162 
163  for(int mem = 0; mem < members ; mem++)
164  {
165  mAnode = mContainer1[mem]/1000;
166  mSeq = mContainer1[mem]%1000;
167 
168  mContainer1[mem] = 0; //zeroed after being used
169 
170  status = hybdata->getListSequences(mAnode,mSequence,sequence);
171 
172  mSeqStart = sequence[mSeq].startTimeBin;
173  mSeqLength = sequence[mSeq].length;
174  mSeqEnd = mSeqStart + mSeqLength - 1;
175 
176  mRightFlagCounter = getSeqOnRight(mAnode,mSeqStart, mSeqEnd, memCount, newMem);
177  mLeftFlagCounter = getSeqOnLeft(mAnode,breakAn,mSeqStart, mSeqEnd, memCount, newMem);
178 
179  // if((mem == members - 1) && (mRightFlagCounter != 0 || mLeftFlagCounter != 0))
180  if((mem == members - 1) && newMem!=0)
181  {
182  members = newMem;
183  mRightFlagCounter = 0;
184  mLeftFlagCounter = 0;
185 
186  // we don't need to do this now since we made each zero already
187 
188  //for(mem = 0; mem < 500; mem++)
189  // mContainer1[mem] = 0;
190 
191  int j = 0;
192 
193  for(mem = 0; mem < members && j<newMem; mem++)
194  {
195  mContainer1[mem] = mContainer2[j];
196  ++j;
197  }
198 
199  for( j = 0; j < newMem; j++)
200  mContainer2[j] = 0;
201 
202  newMem = 0;
203  mem = -1;
204  }
205 
206  // else if((mem == members - 1) && (mRightFlagCounter == 0 && mLeftFlagCounter == 0))
207  else if((mem == members - 1) && newMem==0)
208  {
209  //cout<< memCount <<' '<<"members found in cluster"<<' '<<cluIndex - 1<<endl;
210 
211  mNumOfCluMem[cluIndex-1] = memCount;
212 
213  if( memCount == 1 && mSeqLength==1){
214  cluIndex--;
215  }
216 
217  /* we don't need to do this now since we made each one zero already
218 
219  for(int j = 0; j < 500; j++)
220  {
221  mContainer1[j] = 0;
222  mContainer2[j] = 0;
223  }
224  */
225 
226  mSeq = -1;
227  mAnode = breakAn;
228  status = hybdata->getListSequences(mAnode,mSequence,sequence);
229 
230  }
231  }
232 
233  }
234 
235 //______________________________________________________________________________
236 
237 int StSvtClusterFinder::getSeqOnRight(int mAnode, int mSeqStart, int mSeqEnd, int& memCount,int& newmem)
238 {
239  int mRightAnode = 0, mRightSequence = 0, rightFlagCounter = 0, status;
240  int mRightSeqStart, mRightSeqLength, mRightSeqEnd;
241 
242  mRightAnode = mAnode + 1;
243  if(mRightAnode != 240)
244  {
245  status = hybdata->getListSequences(mRightAnode,mRightSequence,sequence);
246 
247  for(int mRightSeq = 0; mRightSeq < mRightSequence ; mRightSeq++)
248  {
249  if(mAnolist[mAnode] == 240)
250  break;
251 
252  if(mAnolist[mRightAnode] != mAnolist[mAnode] + 1) //are they consecutive?
253  break;
254 
255  mRightSeqStart = sequence[mRightSeq].startTimeBin;
256  mRightSeqLength = sequence[mRightSeq].length;
257  mRightSeqEnd = mRightSeqStart + mRightSeqLength - 1;
258 
259  // cout << "initial clusterflag is " <<' '<<mSeqFlag[mRightAnode][mRightSeq]<<endl;
260 
261  if(mSeqFlag[mRightAnode][mRightSeq] == 0)
262  {
263  if(mRightSeqStart == mSeqStart)
264  {
265  mCluster[cluIndex-1][memCount] = 1000*mRightAnode + mRightSeq;
266  mSeqFlag[mRightAnode][mRightSeq] = 1;
267  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
268  ++memCount;
269  ++newmem;
270  ++rightFlagCounter;
271  }
272 
273  if((mRightSeqStart < mSeqStart)&&(mRightSeqEnd >= mSeqStart - 1))
274  {
275  mCluster[cluIndex-1][memCount] = 1000*mRightAnode + mRightSeq;
276  mSeqFlag[mRightAnode][mRightSeq] = 1;
277  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
278  ++memCount;
279  ++newmem;
280  ++rightFlagCounter;
281  }
282 
283  if((mRightSeqStart > mSeqStart)&&(mRightSeqStart <= mSeqEnd + 1))
284  {
285  mCluster[cluIndex-1][memCount] = 1000*mRightAnode + mRightSeq;
286  mSeqFlag[mRightAnode][mRightSeq] = 1;
287  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
288  ++memCount;
289  ++newmem;
290  ++rightFlagCounter;
291  }
292 
293  //cout<<"mContainer2[newmem-1] ="<<' '<<mContainer2[newmem-1]<<endl;
294  //cout<<"mCluster[cluIndex-1][memCount-1] = "<<' '<<mCluster[cluIndex-1][memCount]<<endl;
295  //cout << "final sequence flag is " <<' '<<mSeqFlag[mRightAnode][mRightSeq]<<endl;
296  //cout << "right flag counter is " <<' '<<rightFlagCounter<<endl;
297  // cout << "new members ="<<' '<< k <<endl;
298  // cout <<"done with mRightSeq"<<endl;
299 
300 
301 
302  }
303  }
304  }
305  return rightFlagCounter;
306 }
307 
308  //_____________________________________________________________________________
309 
310 
311 int StSvtClusterFinder::getSeqOnLeft(int mAnode, int& breakAn,int mSeqStart, int mSeqEnd, int& memCount,int& newmem)
312 {
313  int mLeftAnode = 0, mLeftSequence = 0, leftFlagCounter = 0, status;
314  int mLeftSeqStart, mLeftSeqLength, mLeftSeqEnd;
315 
316  mLeftAnode = mAnode - 1;
317 
318  if(mLeftAnode != -1)
319  {
320  status = hybdata->getListSequences(mLeftAnode,mLeftSequence,sequence);
321 
322  for(int mLeftSeq = 0; mLeftSeq < mLeftSequence ; mLeftSeq++)
323  {
324  if(mAnolist[mAnode] == 1)
325  {
326  breakAn = 0;
327  break;
328  }
329 
330  if(mAnolist[mLeftAnode] != mAnolist[mAnode] - 1) //are they consecutive?
331  {
332  breakAn = mAnode;
333  break;
334  }
335 
336  mLeftSeqStart = sequence[mLeftSeq].startTimeBin;
337  mLeftSeqLength = sequence[mLeftSeq].length;
338  mLeftSeqEnd = mLeftSeqStart + mLeftSeqLength - 1;
339 
340  // cout << "initial clusterflag is " <<' '<<mSeqFlag[mLeftAnode][mLeftSeq]<<endl;
341  if(mSeqFlag[mLeftAnode][mLeftSeq] == 0)
342  {
343  // cout<<"looking for mSeqLeft"<<endl;
344 
345  if(mLeftSeqStart == mSeqStart)
346  {
347  mCluster[cluIndex-1][memCount] = 1000*mLeftAnode + mLeftSeq;
348  mSeqFlag[mLeftAnode][mLeftSeq] = 1;
349  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
350  ++memCount;
351  ++newmem;
352  ++leftFlagCounter;
353  }
354 
355  if((mLeftSeqStart < mSeqStart) && (mLeftSeqEnd >= mSeqStart - 1))
356  {
357  mCluster[cluIndex-1][memCount] = 1000*mLeftAnode + mLeftSeq;
358  mSeqFlag[mLeftAnode][mLeftSeq] = 1;
359  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
360  ++memCount;
361  ++newmem;
362  ++leftFlagCounter;
363  }
364 
365  if((mLeftSeqStart > mSeqStart) && (mLeftSeqStart <= mSeqEnd + 1))
366  {
367  mCluster[cluIndex-1][memCount] = 1000*mLeftAnode + mLeftSeq;
368  mSeqFlag[mLeftAnode][mLeftSeq] = 1;
369  mContainer2[newmem] = mCluster[cluIndex-1][memCount];
370  ++memCount;
371  ++newmem;
372  ++leftFlagCounter;
373  }
374 
375  //cout<<"mContainer2[newmem-1] ="<<' '<<mContainer2[newmem-1]<<endl;
376  //cout<<"mCluster[cluIndex-1][memCount]= "<<' '<<mCluster[cluIndex][memCount]<<endl;
377  //cout << "final sequence flag is " <<' '<<mSeqFlag[mLeftAnode][mLeftSeq]<<endl;
378  //cout << "left flag counter is " <<' '<<leftFlagCounter<<endl;
379  // cout << "done with mLeftSeq"<<endl;
380 
381  }
382  }
383  }
384 
385  return leftFlagCounter;
386 }
387 
388 
389 //____________________________________________________________________________
390 
391 void StSvtClusterFinder::SetHybIndex(int index)
392 {
393  m_hybIndex = index;
394 }
395 
396 int StSvtClusterFinder::ClusterIndex()
397  {
398  return mNumOfClusters;
399  }
400 
401 int StSvtClusterFinder::ClusterMembers(int clu)
402  {
403  return mNumOfCluMem[clu];
404  }
405 
406 int StSvtClusterFinder::ClusterListAnode(int clu, int mem)
407  {
408  return mCluster[clu][mem]/1000;
409  }
410 
411 int StSvtClusterFinder::ClusterSequence(int clu, int mem)
412  {
413  return mCluster[clu][mem]%1000;
414 
415  }
416 
417 int StSvtClusterFinder::ClusterActualAnode(int listAn)
418 {
419  return mAnolist[listAn];
420 }
421 
422 
423 //---------------------------------------------------------------------------
424 
425 void StSvtClusterFinder::ResetContainers()
426 {
427  int i,j;
428  /*
429  for( j = 0; j < 500; j++)
430  {
431  mContainer1[j] = 0;
432  mContainer2[j] = 0;
433  }
434  for ( i=0; i<8000; i++)
435  for ( j=0; j<500; j++)
436  mCluster[i][j] = 0;
437 
438  for ( i=0; i<240; i++)
439  for ( j=0; j<128; j++)
440  mSeqFlag[i][j] = 0;
441  */
442 
443 //array memebers > mNumOfClusters already zero
444 
445  for ( i=0; i< mNumOfClusters; i++)
446  for ( j=0; j< mNumOfCluMem[i]; j++)
447  mCluster[i][j] = 0;
448 
449 // do only the ones that have been flagged
450 
451  {for(int i = 0; i<mNumOfAnodes; i++)
452  for( j=0; j<mNumSeq[i];j++)
453  mSeqFlag[i][j] = 0;}
454 
455  cluIndex=0;
456 }
457 
458 
459 
460