8fe6341f96292a84443969365fdda472afc4aa01
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCPadArray.cxx
1 // @(#) $Id$
2
3 /**************************************************************************
4  * This file is property of and copyright by the ALICE HLT Project        * 
5  * ALICE Experiment at CERN, All rights reserved.                         *
6  *                                                                        *
7  * Authors: Kenneth Aamodt <Kenneth.Aamodt@student.uib.no>                *
8  *          for The ALICE HLT Project.                                    *
9  *                                                                        *
10  * Permission to use, copy, modify and distribute this software and its   *
11  * documentation strictly for non-commercial purposes is hereby granted   *
12  * without fee, provided that the above copyright notice appears in all   *
13  * copies and that both the copyright notice and this permission notice   *
14  * appear in the supporting documentation. The authors make no claims     *
15  * about the suitability of this software for any purpose. It is          *
16  * provided "as is" without express or implied warranty.                  *
17  **************************************************************************/
18
19 /** @file   AliHLTTPCPadArray.cxx
20     @author Kenneth Aamodt
21     @date   
22     @brief  Class containing TPC Pad objects.
23 */
24
25 #if __GNUC__>= 3
26 using namespace std;
27 #endif
28
29 #include <cerrno>
30 #include "AliHLTTPCPadArray.h"
31 #include "AliHLTTPCPad.h"
32 #include "AliHLTStdIncludes.h"
33 #include "AliHLTTPCTransform.h"
34 #include "AliHLTTPCDigitReader.h"
35 #include "AliHLTTPCClusters.h"
36 #include <vector>
37 #include <sys/time.h>
38
39 /** ROOT macro for the implementation of ROOT specific class methods */
40 ClassImp(AliHLTTPCPadArray)
41
42 AliHLTTPCPadArray::AliHLTTPCPadArray()
43   :
44   fPatch(-1),
45   fFirstRow(-1),
46   fLastRow(-1),
47   fThreshold(10),
48   fNumberOfRows(0),
49   fNumberOfPadsInRow(NULL),
50   fDigitReader(NULL),
51   fRowPadVector(0),
52   fClusters(0)
53 {
54   // see header file for class documentation
55   // or
56   // refer to README to build package
57   // or
58   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
59 }
60
61 AliHLTTPCPadArray::AliHLTTPCPadArray(Int_t patch)
62   :
63   fPatch(patch),
64   fFirstRow(-1),
65   fLastRow(-1),
66   fThreshold(10),
67   fNumberOfRows(0),
68   fNumberOfPadsInRow(NULL),
69   fDigitReader(NULL),
70   fRowPadVector(0),
71   fClusters(0)
72 {
73   // see header file for class documentation
74 }
75
76 AliHLTTPCPadArray::AliHLTTPCPadArray(const AliHLTTPCPadArray& srcPadArray)
77   :
78   fPatch(srcPadArray.fPatch),
79   fFirstRow(srcPadArray.fFirstRow),
80   fLastRow(srcPadArray.fLastRow),
81   fThreshold(srcPadArray.fThreshold),
82   fNumberOfRows(srcPadArray.fNumberOfRows),
83   fNumberOfPadsInRow(srcPadArray.fNumberOfPadsInRow),
84   fDigitReader(srcPadArray.fDigitReader),
85   fRowPadVector(srcPadArray.fRowPadVector),
86   fClusters(srcPadArray.fClusters)
87 {
88   // see header file for class documentation
89   HLTFatal("copy constructor not implemented");
90 }
91
92 AliHLTTPCPadArray& AliHLTTPCPadArray::operator=(const AliHLTTPCPadArray&)
93 {
94   // see header file for class documentation
95   HLTFatal("assignment operator not implemented");
96   return (*this);
97 }
98
99 AliHLTTPCPadArray::~AliHLTTPCPadArray()
100 {
101   // see header file for class documentation
102 }
103
104 Int_t AliHLTTPCPadArray::InitializeVector(){
105   // see header file for class documentation
106
107   if(fPatch>5||fPatch<0){
108     HLTFatal("Patch is not set");
109     return 0;
110   }
111
112   fFirstRow = AliHLTTPCTransform::GetFirstRow(fPatch);
113   fLastRow = AliHLTTPCTransform::GetLastRow(fPatch);
114
115   fNumberOfRows=fLastRow-fFirstRow+1;
116   fNumberOfPadsInRow= new Int_t[fNumberOfRows];
117
118   memset( fNumberOfPadsInRow, 0, sizeof(Int_t)*(fNumberOfRows));
119
120   for(Int_t i=0;i<fNumberOfRows;i++){
121     fNumberOfPadsInRow[i]=AliHLTTPCTransform::GetNPads(i+fFirstRow);
122     fPadVector tmpRow;
123     for(Int_t j=0;j<fNumberOfPadsInRow[i];j++){
124       AliHLTTPCPad *tmpPad = new AliHLTTPCPad();
125       tmpPad->SetID(i,j);
126       tmpRow.push_back(tmpPad);
127     }
128     fRowPadVector.push_back(tmpRow);
129   }
130 }
131
132 Int_t AliHLTTPCPadArray::DeInitializeVector(){
133   for(Int_t i=0;i<fNumberOfRows;i++){
134     for(Int_t j=0;j<fNumberOfPadsInRow[i];j++){
135       delete fRowPadVector[i][j];
136     }
137     fRowPadVector[i].clear();
138   }
139   fRowPadVector.clear();
140   return 1;
141
142 void AliHLTTPCPadArray::SetPatch(Int_t patch){
143   fPatch=patch;
144 }
145
146 void AliHLTTPCPadArray::SetDigitReader(AliHLTTPCDigitReader* digitReader){
147   fDigitReader=digitReader;
148 }
149 Int_t AliHLTTPCPadArray::ReadData(){
150
151   switch (fPatch){
152   case 0:
153     while(fDigitReader->Next()){
154       fRowPadVector[fDigitReader->GetRow()-fFirstRow][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());
155     }
156     break;
157   case 1:
158     while(fDigitReader->Next()){
159       fRowPadVector[fDigitReader->GetRow()-fFirstRow][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());      
160     }
161     break;
162   case 2:
163     while(fDigitReader->Next()){
164       fRowPadVector[fDigitReader->GetRow()][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());
165     }
166     break;
167   case 3:
168     while(fDigitReader->Next()){
169       fRowPadVector[fDigitReader->GetRow()-27][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());
170     }
171     break;
172   case 4:
173     while(fDigitReader->Next()){
174       fRowPadVector[fDigitReader->GetRow()-54][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());
175     }
176     break;
177   case 5:
178     while(fDigitReader->Next()){
179       fRowPadVector[fDigitReader->GetRow()-76][fDigitReader->GetPad()]->SetDataSignal(fDigitReader->GetTime(),fDigitReader->GetSignal());
180     }
181     break;
182   }
183   return 0;
184 }
185 void AliHLTTPCPadArray::FindClusters(Int_t match){
186   //see header file for documentation
187   Int_t nClusters=0;
188   Int_t totalChargeOfPreviousCandidate=0;
189   Int_t clusterChargeIsFalling=0;
190   for(Int_t row=0;row<fNumberOfRows;row++){
191     for(Int_t pad=0;pad<fNumberOfPadsInRow[row]-1;pad++){
192       AliHLTTPCPad *tmp1=fRowPadVector[row][pad];
193       AliHLTTPCPad *tmp2=fRowPadVector[row][pad+1];
194       for(Int_t c1=0;c1<tmp1->fClusterCandidates.size();c1++){
195
196         if(tmp1->fUsedClusterCandidates[c1]){
197           continue;
198         }
199
200         for(Int_t c2=0;c2<tmp2->fClusterCandidates.size();c2++){
201
202           if(tmp2->fUsedClusterCandidates[c2]){
203             continue;
204           }
205
206           Int_t diff= tmp1->fClusterCandidates[c1].fMean - tmp2->fClusterCandidates[c2].fMean;
207
208           if(diff < -match){
209             break;
210           }
211           if(diff <= match){
212
213             if((Int_t)(tmp1->fClusterCandidates[c1].fTotalCharge - tmp2->fClusterCandidates[c2].fTotalCharge)>0){
214               clusterChargeIsFalling=1;
215             }
216
217             tmp1->fUsedClusterCandidates[c1]=1;
218             tmp2->fUsedClusterCandidates[c2]=1;
219
220             AliHLTTPCClusters tmpCluster;
221             tmpCluster.fTotalCharge = tmp1->fClusterCandidates[c1].fTotalCharge;
222             tmpCluster.fPad     = tmp1->fClusterCandidates[c1].fPad;
223             tmpCluster.fPad2    = tmp1->fClusterCandidates[c1].fPad2;
224             tmpCluster.fTime    = tmp1->fClusterCandidates[c1].fTime;
225             tmpCluster.fTime2   = tmp1->fClusterCandidates[c1].fTime2;
226             tmpCluster.fRowNumber = row;
227
228             tmpCluster.fTotalCharge += tmp2->fClusterCandidates[c2].fTotalCharge;
229             tmpCluster.fPad     += tmp2->fClusterCandidates[c2].fPad;
230             tmpCluster.fPad2    += tmp2->fClusterCandidates[c2].fPad2;
231             tmpCluster.fTime    += tmp2->fClusterCandidates[c2].fTime;
232             tmpCluster.fTime2   += tmp2->fClusterCandidates[c2].fTime2;
233             tmpCluster.fMean         = tmp2->fClusterCandidates[c2].fMean;
234             totalChargeOfPreviousCandidate = tmp2->fClusterCandidates[c2].fTotalCharge;
235
236             int rowNumber=row;
237             int lastPad=pad+1;
238             nClusters++;
239             Int_t doBreak=0;
240             for(Int_t morePads=pad+2;morePads<fNumberOfPadsInRow[row];morePads++){
241               AliHLTTPCPad *tmpx=fRowPadVector[row][morePads];
242               if(morePads>lastPad+1){
243                 break;
244               }
245               for(Int_t cx=0;cx<tmpx->fClusterCandidates.size();cx++){
246                 if(tmpx->fUsedClusterCandidates[cx]){
247                   continue;
248                 }
249                 Int_t diffx=tmpCluster.fMean - tmpx->fClusterCandidates[cx].fMean;
250                 if(diffx<-match){
251                   doBreak=1;
252                   break;
253                 }
254                 if(diffx <= match){
255                   if((Int_t)(totalChargeOfPreviousCandidate - tmpx->fClusterCandidates[cx].fTotalCharge)>0){
256                     clusterChargeIsFalling=1;
257                   }
258
259                   if(clusterChargeIsFalling&&(Int_t)(totalChargeOfPreviousCandidate - tmpx->fClusterCandidates[cx].fTotalCharge)<=0){
260                     //Means we have a deconvoluted cluster.
261                     totalChargeOfPreviousCandidate=0;
262                     doBreak=1;
263                     break;
264                   }
265                   
266                   tmpx->fUsedClusterCandidates[cx]=1;
267                   tmpCluster.fTotalCharge += tmpx->fClusterCandidates[cx].fTotalCharge;
268                   tmpCluster.fPad     += tmpx->fClusterCandidates[cx].fPad;
269                   tmpCluster.fPad2    += tmpx->fClusterCandidates[cx].fPad2;
270                   tmpCluster.fTime    += tmpx->fClusterCandidates[cx].fTime;
271                   tmpCluster.fTime2   += tmpx->fClusterCandidates[cx].fTime2;
272                   tmpCluster.fMean         = tmpx->fClusterCandidates[cx].fMean;
273                   lastPad=morePads;
274
275                   totalChargeOfPreviousCandidate=tmpx->fClusterCandidates[cx].fTotalCharge;
276                 }
277               }
278               if(doBreak){
279                 break;
280               }
281             }
282             
283             if(tmpCluster.fTotalCharge<fThreshold){
284               nClusters--;
285             }
286             else{
287               //Code to look for tails, TODO insert flag.
288               /* UInt_t meanTime=tmpCluster.fMean;
289               if(pad>0){
290                 AliHLTTPCPad *tmpBefore=fRowPadVector[row][pad-1];
291                 //checking the fMean -1 timebin for single timebin value in the pad before the cluster
292                 if(meanTime>0){
293                   Int_t charge =tmpBefore->GetDataSignal(meanTime-1); 
294                   if(charge){
295                     tmpCluster.fTotalCharge+= charge;
296                     tmpCluster.fPad        += charge*(pad-1);
297                     tmpCluster.fPad2       += charge*(pad-1)*(pad-1); 
298                     tmpCluster.fTime       += meanTime*charge;
299                     tmpCluster.fTime2      += meanTime*charge*charge;
300                   } 
301                 }
302                 //checking the fMean timebin for single timebin value in the pad before the cluster
303                 Int_t charge2 =tmpBefore->GetDataSignal(meanTime); 
304                 if(charge2){
305                   tmpCluster.fTotalCharge+= charge2;
306                   tmpCluster.fPad        += charge2*(pad);
307                   tmpCluster.fPad2       += charge2*(pad)*(pad); 
308                   tmpCluster.fTime       += meanTime*charge2;
309                   tmpCluster.fTime2      += meanTime*charge2*charge2;
310                 } 
311                 //checking the fMean +1 timebin for single timebin value in the pad before the cluster
312                 if(meanTime<AliHLTTPCTransform::GetNTimeBins()){
313                   Int_t charge3 =tmpBefore->GetDataSignal(meanTime+1); 
314                   if(charge3){
315                     tmpCluster.fTotalCharge+= charge3;
316                     tmpCluster.fPad        += charge3*(pad+1);
317                     tmpCluster.fPad2       += charge3*(pad+1)*(pad+1); 
318                     tmpCluster.fTime       += meanTime*charge3;
319                     tmpCluster.fTime2      += meanTime*charge3*charge3;
320                   } 
321                 }
322               }
323               
324               if(lastPad<fNumberOfPadsInRow[row]-2){
325                 AliHLTTPCPad *tmpAfter=fRowPadVector[row][lastPad+1];
326                 //checking the fMean +1 timebin for single timebin value in the pad after the cluster
327                 if(meanTime>0){
328                   Int_t charge4 =tmpAfter->GetDataSignal(meanTime-1); 
329                   if(charge4){
330                     tmpCluster.fTotalCharge+= charge4;
331                     tmpCluster.fPad        += charge4*(pad-1);
332                     tmpCluster.fPad2       += charge4*(pad-1)*(pad-1); 
333                     tmpCluster.fTime       += meanTime*charge4;
334                     tmpCluster.fTime2      += meanTime*charge4*charge4;
335                   } 
336                 }
337
338               
339                 //checking the fMean +1 timebin for single timebin value in the pad after the cluster
340                 Int_t charge5 =tmpAfter->GetDataSignal(meanTime); 
341                 if(charge5){
342                   tmpCluster.fTotalCharge+= charge5;
343                   tmpCluster.fPad        += charge5*(pad);
344                   tmpCluster.fPad2       += charge5*(pad)*(pad); 
345                   tmpCluster.fTime       += meanTime*charge5;
346                   tmpCluster.fTime2      += meanTime*charge5*charge5;
347                 } 
348                 //checking the fMean +1 timebin for single timebin value in the pad after the cluster
349                 if(meanTime<AliHLTTPCTransform::GetNTimeBins()){
350                   Int_t charge6 =tmpAfter->GetDataSignal(meanTime+1); 
351                   if(charge6){
352                     tmpCluster.fTotalCharge+= charge6;
353                     tmpCluster.fPad        += charge6*(pad+1);
354                     tmpCluster.fPad2       += charge6*(pad+1)*(pad+1); 
355                     tmpCluster.fTime       += meanTime*charge6;
356                     tmpCluster.fTime2      += meanTime*charge6*charge6;
357                   } 
358                 }
359               }
360               */
361               //              tmpCluster.fTime= tmpCluster.fTime/tmpCluster.fTotalCharge;
362               totalChargeOfPreviousCandidate=0;
363               clusterChargeIsFalling=0;
364               tmpCluster.fFirstPad=pad;
365               switch (fPatch){
366               case 0:
367                 tmpCluster.fRowNumber=row;
368                 break;
369               case 1:
370                 tmpCluster.fRowNumber=row+30;
371                 break;
372               case 2:
373                 tmpCluster.fRowNumber=row+63;
374                 break;
375               case 3:
376                 tmpCluster.fRowNumber=row+90;
377                 break;
378               case 4:
379                 tmpCluster.fRowNumber=row+117;
380                 break;
381               case 5:
382                 tmpCluster.fRowNumber=row+139;
383                 break;
384               }
385
386               fClusters.push_back(tmpCluster);
387             }
388           }
389         }
390       }
391     }
392   }
393
394   HLTInfo("Found %d clusters.",nClusters);
395   // PrintClusters();
396 }
397 void AliHLTTPCPadArray::PrintClusters(){
398   for(int i=0;i<fClusters.size();i++){
399     cout<<"Cluster number: "<<i<<endl;
400     cout<<"Row: "<<fClusters[i].fRowNumber<<" Pad: "<<fClusters[i].fFirstPad<<endl;
401     cout<<"Total Charge: "<<fClusters[i].fTotalCharge<<endl;
402     cout<<"PadError:     "<<fClusters[i].fPad2<<endl;
403     cout<<"TimeMean:     "<<fClusters[i].fTime<<endl;
404     cout<<"TimeError:    "<<fClusters[i].fTime2<<endl;
405     cout<<endl;
406     cout<<endl;
407   }
408 }