Updated version from B. Batyunya and E. Fragiacomo
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderSSD.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 #include <iostream.h>
17 #include <TArrayI.h>
18 #include "AliRun.h"
19 #include "AliITS.h"
20 #include "AliITSdigit.h"
21 #include "AliITSRawCluster.h"
22 #include "AliITSRecPoint.h"
23 #include "AliITSMapA1.h"
24 #include "AliITSClusterFinderSSD.h"
25 #include "AliITSclusterSSD.h"
26 #include "AliITSpackageSSD.h"
27 #include "AliITSsegmentation.h"
28 #include "AliITSgeom.h"
29
30 const Bool_t AliITSClusterFinderSSD::fgkSIDEP=kTRUE;
31 const Bool_t AliITSClusterFinderSSD::fgkSIDEN=kFALSE;
32 const Int_t debug=0;
33
34 ClassImp(AliITSClusterFinderSSD)
35
36 //____________________________________________________________________
37 //
38 //  Constructor
39 //____________________________________________________________________
40 //                                     
41
42
43 AliITSClusterFinderSSD::AliITSClusterFinderSSD(AliITSsegmentation *seg, TClonesArray *digits)   
44 {
45 //Standard constructor
46
47     fSegmentation=seg;
48     fDigits=digits;
49     
50     fMap = new AliITSMapA1(fSegmentation,fDigits);  
51   
52     fITS=(AliITS*)gAlice->GetModule("ITS");
53     
54     fClusterP =  new TClonesArray ("AliITSclusterSSD",200);    
55     fNClusterP =0;   
56     
57     fClusterN=  new TClonesArray ("AliITSclusterSSD",200);   
58     fNClusterN =0; 
59
60     fPackages =  new TClonesArray ("AliITSpackageSSD",200);    //packages  
61     fNPackages = 0;
62
63         
64     fDigitsIndexP      =  new TArrayI(300);
65     fNDigitsP          =  0;
66     
67     fDigitsIndexN      =  new TArrayI(300);
68     fNDigitsN          =  0;
69     
70     fPitch = fSegmentation->Dpx(0);
71     fPNsignalRatio=7./8.;         // warning: hard-wired number
72
73 }
74
75 //-------------------------------------------------------
76 AliITSClusterFinderSSD::~AliITSClusterFinderSSD() {
77 // Default destructor   
78
79     delete fClusterP;
80     delete fClusterN;        
81     delete fPackages;        
82     delete fDigitsIndexP;        
83     delete fDigitsIndexN; 
84     delete fMap;
85
86 }
87
88 //-------------------------------------------------------
89 void AliITSClusterFinderSSD::InitReconstruction()
90 {
91 // initialization of the cluster finder
92
93   register Int_t i; //iterator
94
95   for (i=0;i<fNClusterP;i++)
96     {
97       fClusterP->RemoveAt(i);
98     }
99   fNClusterP  =0;
100   for (i=0;i<fNClusterN;i++)
101     {
102       fClusterN->RemoveAt(i);
103     }
104   fNClusterN=0;
105
106   for (i=0;i<fNPackages;i++)
107     {
108       fPackages->RemoveAt(i);
109     }
110
111   fNPackages = 0;
112   fNDigitsP=0;
113   fNDigitsN=0;
114
115   Float_t StereoP,StereoN;
116   fSegmentation->Angles(StereoP,StereoN);
117
118   CalcStepFactor (StereoP,StereoN);
119
120   if (debug) cout<<"fSFF = "<<fSFF<<"  fSFB = "<<fSFB<<"\n";
121 }
122
123
124 //---------------------------------------------
125 void AliITSClusterFinderSSD::FindRawClusters(Int_t module) 
126 {
127 // This function findes out all clusters belonging to one module
128 // 1. Zeroes all space after previous module reconstruction
129 // 2. Finds all neighbouring digits, create clusters
130 // 3. If necesery, resolves for each group of neighbouring digits 
131 //    how many clusters creates it.
132 // 4. Colculate the x and z coordinate  
133
134   Int_t lay, lad, detect;
135   AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
136   AliITSgeom *geom = aliITS->GetITSgeom();
137             geom->GetModuleId(module,lay, lad, detect);
138             if ( lay == 6 )
139               ((AliITSsegmentationSSD*)fSegmentation)->SetLayer(6);
140             if ( lay == 5 )
141               ((AliITSsegmentationSSD*)fSegmentation)->SetLayer(5);
142
143   InitReconstruction();  //ad. 1
144   fMap->FillMap();
145   FillDigitsIndex();
146   SortDigits();
147   FindNeighbouringDigits(); //ad. 2
148   //SeparateOverlappedClusters();  //ad. 3
149   ClustersToPackages();  //ad. 4
150   AliITSRecPoint rnew;
151   fMap->ClearMap();
152 }
153
154
155 //-------------------------------------------------
156 void AliITSClusterFinderSSD::FindNeighbouringDigits()
157 {
158 //If there are any digits on this side, create 1st Cluster,
159 // add to it this digit, and increment number of clusters
160
161
162  register Int_t i;
163
164
165  if ((fNDigitsP>0 )  && (fNDigitsN > 0 )) {     
166
167    Int_t currentstripNo;
168    Int_t *dbuffer = new Int_t [300];   //buffer for strip numbers
169    Int_t dnumber;    //curent number of digits in buffer
170    TArrayI      &lDigitsIndexP = *fDigitsIndexP;
171    TArrayI      &lDigitsIndexN = *fDigitsIndexN;
172    TObjArray    &lDigits       = *(Digits());
173    TClonesArray &lClusterP     = *fClusterP;
174    TClonesArray &lClusterN     = *fClusterN;
175   
176    //process P side 
177    dnumber = 1;
178    dbuffer[0]=lDigitsIndexP[0];
179    //If next digit is a neighbour of previous, adds to last cluster this digit
180    for (i=1; i<fNDigitsP; i++) {
181      //reads new digit
182      currentstripNo = ((AliITSdigitSSD*)lDigits[lDigitsIndexP[i]])->
183                                                             GetStripNumber(); 
184      //check if it is a neighbour of a previous one
185      if ( (((AliITSdigitSSD*)lDigits[lDigitsIndexP[i-1]])->GetStripNumber()) 
186            ==  (currentstripNo-1) ) dbuffer[dnumber++]=lDigitsIndexP[i];
187      else  { 
188        //create a new one side cluster
189        new(lClusterP[fNClusterP++]) AliITSclusterSSD(dnumber,dbuffer,Digits(),fgkSIDEP); 
190        dbuffer[0]=lDigitsIndexP[i];
191        dnumber = 1;
192      }
193    } // end loop over fNDigitsP
194    new(lClusterP[fNClusterP++]) AliITSclusterSSD(dnumber,dbuffer,Digits(),fgkSIDEP);
195
196
197    //process N side 
198    //for comments, see above
199    dnumber = 1;
200    dbuffer[0]=lDigitsIndexN[0];
201    //If next digit is a neighbour of previous, adds to last cluster this digit
202    for (i=1; i<fNDigitsN; i++) { 
203      currentstripNo = ((AliITSdigitSSD*)(lDigits[lDigitsIndexN[i]]))->
204                                                             GetStripNumber();
205      if ( (((AliITSdigitSSD*)lDigits[lDigitsIndexN[i-1]])->GetStripNumber()) 
206            == (currentstripNo-1) ) dbuffer[dnumber++]=lDigitsIndexN[i];
207      else {
208        new(lClusterN[fNClusterN++]) AliITSclusterSSD(dnumber,dbuffer,Digits(),fgkSIDEN);
209        dbuffer[0]=lDigitsIndexN[i];
210        dnumber = 1;
211      }
212    } // end loop over fNDigitsN
213    new(lClusterN[fNClusterN++]) AliITSclusterSSD(dnumber,dbuffer,Digits(),fgkSIDEN);
214    delete [] dbuffer; 
215  
216  } // end condition on  NDigits 
217
218  if (debug) cout<<"\n Found clusters: fNClusterP = "<<fNClusterP<<"  fNClusterN ="<<fNClusterN<<"\n";
219
220 }  
221
222 //--------------------------------------------------------------
223
224 void AliITSClusterFinderSSD::SeparateOverlappedClusters()
225 {
226 // overlapped clusters separation
227
228   register Int_t i; //iterator
229
230   Float_t  factor=0.75;            // How many percent must be lower signal 
231                                    // on the middle one digit
232                                    // from its neighbours
233   Int_t    signal0;                //signal on the strip before the current one
234   Int_t    signal1;                //signal on the current one signal
235   Int_t    signal2;                //signal on the strip after the current one
236   TArrayI *splitlist;              //  List of splits
237   Int_t    numerofsplits=0;        // number of splits
238   Int_t    initPsize = fNClusterP; //initial size of the arrays 
239   Int_t    initNsize = fNClusterN; //we have to keep it because it will grow 
240                                    // in this function and it doasn't make 
241                                    // sense to pass through it again
242
243   splitlist = new TArrayI(300);
244
245   for (i=0;i<initPsize;i++)
246   {
247     if (( ((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits())==1) continue;
248     if (( ((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits())==2) continue;
249         Int_t nj=(((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits()-1);
250         for (Int_t j=1; j<nj; j++)
251           {
252             signal1=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j);
253             signal0=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j-1);
254             signal2=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j+1);
255             //if signal is less then factor*signal of its neighbours
256             if (  (signal1<(factor*signal0)) && (signal1<(factor*signal2)) )
257              {                                                               
258                (*splitlist)[numerofsplits++]=j;  
259              }
260           } // end loop over number of digits
261           //split this cluster if necessary
262           if(numerofsplits>0) SplitCluster(splitlist,numerofsplits,i,fgkSIDEP);
263           numerofsplits=0;
264
265           //in signed places (splitlist)
266   } // end loop over clusters on Pside
267
268   for (i=0;i<initNsize;i++) {
269     if (( ((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits())==1) continue;
270     if (( ((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits())==2) continue;
271        Int_t nj=(((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits()-1);
272        for (Int_t j=1; j<nj; j++)
273           {
274             signal1=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j);
275             signal0=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j-1);
276             signal2=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j+1);
277             //if signal is less then factor*signal of its neighbours
278             if (  (signal1<(factor*signal0)) && (signal1<(factor*signal2)) ) 
279                (*splitlist)[numerofsplits++]=j;  
280           } // end loop over number of digits 
281           //split this cluster into more clusters
282           if(numerofsplits>0) SplitCluster(splitlist,numerofsplits,i,fgkSIDEN);
283           numerofsplits=0;                                                              //in signed places (splitlist)
284   } // end loop over clusters on Nside
285
286   delete splitlist;
287 }
288
289 //-------------------------------------------------------
290 void AliITSClusterFinderSSD::SplitCluster(TArrayI *list, Int_t nsplits, Int_t index, Bool_t side)
291 {
292   //This function splits one side cluster into more clusters
293   //number of splits is defined by "nsplits"
294   //Place of splits are defined in the TArray "list" 
295   
296   // For further optimisation: Replace this function by two 
297   // specialised ones (each for one side)
298   // save one "if"
299
300   //For comlete comments see AliITSclusterSSD::SplitCluster
301
302
303   register Int_t i; //iterator
304
305   AliITSclusterSSD* curentcluster;
306   Int_t   *tmpdigits = new Int_t[100];
307   Int_t    NN;
308   // side true means P side
309   if (side) {
310      curentcluster =((AliITSclusterSSD*)((*fClusterP)[index])) ;
311      for (i = nsplits; i>0 ;i--) {  
312          NN=curentcluster->SplitCluster((*list)[(i-1)],tmpdigits);
313          new ((*fClusterP)[fNClusterP]) AliITSclusterSSD(NN,tmpdigits,Digits(),side);
314          ( (AliITSclusterSSD*)((*fClusterP)[fNClusterP]) )->
315                                                       SetLeftNeighbour(kTRUE);
316          //if left cluster had neighbour on the right before split 
317          //new should have it too
318          if ( curentcluster->GetRightNeighbour() ) 
319                      ( (AliITSclusterSSD*)((*fClusterP)[fNClusterP]) )->
320                                                      SetRightNeighbour(kTRUE);
321          else curentcluster->SetRightNeighbour(kTRUE); 
322          fNClusterP++;
323      } // end loop over nplits
324   } else {
325      curentcluster =((AliITSclusterSSD*)((*fClusterN)[index]));
326      for (i = nsplits; i>0 ;i--) {  
327          NN=curentcluster->SplitCluster((*list)[(i-1)],tmpdigits);
328          new ((*fClusterN)[fNClusterN]) AliITSclusterSSD(NN,tmpdigits,Digits(),side);
329          ((AliITSclusterSSD*)((*fClusterN)[fNClusterN]))->
330                                                     SetRightNeighbour(kTRUE);
331          if (curentcluster->GetRightNeighbour())
332                       ( (AliITSclusterSSD*)( (*fClusterN)[fNClusterN]) )->
333                                                      SetRightNeighbour(kTRUE);
334          else curentcluster->SetRightNeighbour(kTRUE);      
335          fNClusterN++;
336      } // end loop over nplits
337   } // end if side
338   delete []tmpdigits;
339
340 }
341
342
343 //-------------------------------------------------
344 Int_t AliITSClusterFinderSSD::SortDigitsP(Int_t start, Int_t end)
345 {
346 // sort digits on the P side
347
348   
349   Int_t right;
350   Int_t left;
351   if (start != (end - 1) ) 
352     {
353       left=this->SortDigitsP(start,(start+end)/2);
354       right=this->SortDigitsP((start+end)/2,end);  
355       return (left || right);
356     }
357   else
358    { 
359     left =  ((AliITSdigitSSD*)((*(Digits()))[(*fDigitsIndexP)[start]]))->GetStripNumber();
360     right= ((AliITSdigitSSD*)((*(Digits()))[(*fDigitsIndexP)[end]]))->GetStripNumber();  
361     if( left > right )
362      {
363        Int_t tmp = (*fDigitsIndexP)[start];
364        (*fDigitsIndexP)[start]=(*fDigitsIndexP)[end];
365        (*fDigitsIndexP)[end]=tmp;
366        return 1;
367      }
368     else return 0;
369    }
370 }
371
372
373 //--------------------------------------------------
374
375 Int_t AliITSClusterFinderSSD::SortDigitsN(Int_t start, Int_t end)
376 {
377 // sort digits on the N side
378
379   Int_t right;
380   Int_t left;
381   if (start != (end - 1)) 
382     {
383       left=this->SortDigitsN(start,(start+end)/2);
384       right=this->SortDigitsN((start+end)/2,end);  
385       return (left || right);
386     }
387   else 
388    {
389     left =((AliITSdigitSSD*)((*(Digits()))[(*fDigitsIndexN)[start]]))->GetStripNumber(); 
390     right=((AliITSdigitSSD*)((*(Digits()))[(*fDigitsIndexN)[end]]))->GetStripNumber();
391     if ( left > right )
392       {
393         Int_t tmp = (*fDigitsIndexN)[start];
394         (*fDigitsIndexN)[start]=(*fDigitsIndexN)[end];
395         (*fDigitsIndexN)[end]=tmp;
396         return 1;
397       }else return 0;
398    }  
399 }
400
401
402 //------------------------------------------------
403 void AliITSClusterFinderSSD::FillDigitsIndex()
404 {
405  //Fill the indexes of the clusters belonging to a given ITS module
406
407  Int_t PNs=0, NNs=0;
408  Int_t tmp,bit,k;
409  Int_t N;
410  Int_t i;
411  
412  N = fDigits->GetEntriesFast();
413
414  Int_t* PSidx = new Int_t [N*sizeof(Int_t)];
415  Int_t* NSidx = new Int_t [N*sizeof(Int_t)]; 
416  if (fDigitsIndexP==NULL) fDigitsIndexP = new TArrayI(N);
417  if (fDigitsIndexN==NULL) fDigitsIndexN = new TArrayI(N);
418  
419  AliITSdigitSSD *dig; 
420  
421  for ( i = 0 ; i< N; i++ ) {
422       dig = (AliITSdigitSSD*)GetDigit(i);
423       if(dig->IsSideP()) { 
424            bit=1;
425            tmp=dig->GetStripNumber();
426            // I find this totally unnecessary - it's just a 
427            // CPU consuming double check
428            for( k=0;k<PNs;k++)
429             {
430              if (tmp==PSidx[k])
431               {
432                  if (debug) cout<<"Such a digit exists \n";
433                  bit=0;
434               }
435            }   
436            // end comment 
437         if(bit) {
438             fDigitsIndexP->AddAt(i,fNDigitsP++);
439             PSidx[PNs++]=tmp;
440         }
441       } else {
442          bit=1;
443          tmp=dig->GetStripNumber();
444          // same as above
445          for( k=0;k<NNs;k++)
446           {
447            if (tmp==NSidx[k])
448             {
449              if (debug) cout<<"Such a digit exists \n";
450              bit=0;
451             }
452           } 
453          // end comment
454          if (bit) {
455           fDigitsIndexN->AddAt(i,fNDigitsN++);
456           NSidx[NNs++] =tmp;
457          }
458       }
459  }
460    
461  delete [] PSidx;
462  delete [] NSidx;
463
464  if (debug) cout<<"Digits :  P = "<<fNDigitsP<<"   N = "<<fNDigitsN<<endl;
465
466 }
467
468
469 //-------------------------------------------
470
471 void AliITSClusterFinderSSD::SortDigits()
472 {
473 // sort digits
474
475   Int_t i;
476   if(fNDigitsP>1) 
477   for (i=0;i<fNDigitsP-1;i++)
478   if (SortDigitsP(0,(fNDigitsP-1-i))==0) break;
479
480   if(fNDigitsN>1) 
481     for (i=0;i<fNDigitsN-1;i++)
482   if(SortDigitsN(0,(fNDigitsN-1-i))==0) break;
483 }
484
485
486
487 //----------------------------------------------
488 void AliITSClusterFinderSSD::FillClIndexArrays(Int_t* arrayP, Int_t *arrayN)
489 {
490 // fill cluster index array
491
492   register Int_t i;
493   for (i=0; i<fNClusterP;i++)
494     {
495       arrayP[i]=i;
496     }
497   for (i=0; i<fNClusterN;i++)
498     {
499       arrayN[i]=i;
500     }
501 }
502
503
504 //------------------------------------------------------
505 void AliITSClusterFinderSSD::SortClusters(Int_t* arrayP, Int_t *arrayN)
506 {
507 // sort clusters
508
509   Int_t i;
510   if(fNClusterP>1) 
511     for (i=0;i<fNClusterP-1;i++)
512       if (SortClustersP(0,(fNClusterP-1),arrayP)==0)  break;
513     
514     
515   if(fNClusterN>1) 
516     for (i=0;i<fNClusterN-1;i++)
517       if (SortClustersN(0,(fNClusterN-1),arrayN)==0)  break;
518
519 }
520
521
522
523 //---------------------------------------------------
524 Int_t AliITSClusterFinderSSD::SortClustersP(Int_t start, Int_t end, Int_t *array)
525 {
526 //Sort P side clusters
527
528
529   Int_t right;
530   Int_t left;
531   if (start != (end - 1) ) {
532       left=this->SortClustersP(start,(start+end)/2,array);
533       right=this->SortClustersP((start+end)/2,end,array);  
534       return (left || right);
535   } else {
536       left =((AliITSclusterSSD*)((*fClusterP)[array[start]]))->
537                                                          GetDigitStripNo(0);
538       right=((AliITSclusterSSD*)((*fClusterP)[array[ end ]]))->
539                                                          GetDigitStripNo(0);
540       if(left>right) {
541          Int_t tmp = array[start];
542          array[start]=array[end];
543          array[end]=tmp;
544          return 1;
545       } else return 0;
546   }
547
548  
549 }
550
551
552
553 //-------------------------------------------------------
554 Int_t AliITSClusterFinderSSD::SortClustersN(Int_t start, Int_t end, Int_t *array)
555 {
556 //Sort N side clusters
557
558
559   Int_t right;
560   Int_t left;
561   
562   if (start != (end - 1) ) {
563       left=this->SortClustersN(start,(start+end)/2,array);
564       right=this->SortClustersN((start+end)/2,end,array);  
565       return (left || right);
566   } else {
567       left =((AliITSclusterSSD*)((*fClusterN)[array[start]]))->
568                                                          GetDigitStripNo(0);
569       right=((AliITSclusterSSD*)((*fClusterN)[array[ end ]]))->
570                                                          GetDigitStripNo(0);
571       if( left > right) {
572          Int_t tmp = array[start];
573          array[start]=array[end];
574          array[end]=tmp;
575          return 1;
576       } else return 0;
577     } 
578
579 }
580     
581
582
583 //-------------------------------------------------------
584 void AliITSClusterFinderSSD::ClustersToPackages()
585 {  
586 // fill packages   
587     
588   Int_t *oneSclP = new Int_t[fNClusterP]; //I want to have sorted 1 S clusters
589   Int_t *oneSclN = new Int_t[fNClusterN]; //I can not sort it in TClonesArray
590                                           //so, I create table of indexes and 
591                                           //sort it
592                                           //I do not use TArrayI on purpose
593                                           // MB: well, that's not true that one
594                                           //cannot sort objects in TClonesArray
595   
596   //AliITSpackageSSD *currentpkg;
597   AliITSclusterSSD *currentP;
598   AliITSclusterSSD *currentN;
599   Int_t j1, j2;    
600
601 //Fills in One Side Clusters Index Array
602   FillClIndexArrays(oneSclP,oneSclN); 
603 //Sorts filled Arrays
604   SortClusters(oneSclP,oneSclN);                   
605
606   fNPackages=1;      
607   new ((*fPackages)[0]) AliITSpackageSSD(fClusterP,fClusterN);
608   //currentpkg = (AliITSpackageSSD*)((*fPackages)[0]);
609
610   // Take all pairs of recpoints x coordinates in both sides  
611   // to calculate z coordinates of the recpoints
612   for (j1=0;j1<fNClusterP;j1++) {  
613       currentP = GetPSideCluster(oneSclP[j1]);
614       Double_t xP = currentP->GetPosition();
615       //Int_t NxP = currentP->GetNumOfDigits();
616       Float_t signalP = currentP->GetTotalSignal();
617     for (j2=0;j2<fNClusterN;j2++) {  
618       currentN = GetNSideCluster(oneSclN[j2]);
619       Double_t xN = currentN->GetPosition();
620       //Int_t NxN = currentN->GetNumOfDigits();
621       Float_t signalN = currentN->GetTotalSignal();
622       CreateNewRecPoint(xP,1,xN,1,signalP,signalN,currentP, currentN, 0.75);
623
624     }
625   }
626
627    delete oneSclP;
628    delete oneSclN;
629
630 }
631
632
633 //------------------------------------------------------
634 Bool_t AliITSClusterFinderSSD::    
635 CreateNewRecPoint(Float_t P, Float_t dP, Float_t N, Float_t dN,
636                   Float_t SigP,Float_t SigN, 
637                   AliITSclusterSSD *clusterP, AliITSclusterSSD *clusterN,
638                   Stat_t prob)
639 {
640 // create the recpoints
641   const Float_t kADCtoKeV = 2.16; 
642   // 50 ADC units -> 30000 e-h pairs; 1e-h pair -> 3.6e-3 KeV;
643   // 1 ADC unit -> (30000/50)*3.6e-3 = 2.16 KeV 
644   const Float_t kconv = 1.0e-4; 
645
646   const Float_t kRMSx = 20.0*kconv; 
647   const Float_t kRMSz = 800.0*kconv; 
648
649   Int_t n=0;
650   Int_t *tr;
651   Int_t NTracks;
652
653   if (GetCrossing(P,N)) {
654   
655     //GetCrossingError(dP,dN);
656      AliITSRawClusterSSD cnew;
657      Int_t nstripsP=clusterP->GetNumOfDigits();
658      Int_t nstripsN=clusterN->GetNumOfDigits();
659      Float_t signal = 0;
660      Float_t dedx = 0;
661
662      if(SigP>SigN) {
663        signal = SigP;
664        dedx = SigP*kADCtoKeV;
665      }else{
666        signal = SigN;
667        dedx = SigN*kADCtoKeV;
668      }         
669
670      Float_t signalCut = TMath::Abs((SigP-SigN)/signal);
671
672
673      tr = (Int_t*) clusterP->GetTracks(n);
674      NTracks = clusterP->GetNTracks();
675      //cout<<"NewRec: NTracks,tr0,tr1,tr2 ="<<NTracks<<","<<tr[0]<<","<<tr[1]<<","<<tr[2]<<endl;
676      if(NTracks>0 && signalCut<0.18) {
677      // the cut of the signals difference in P and N sides to
678      // remove the "ghosts"
679        cnew.fSignalP=SigP;
680        cnew.fSignalN=SigN;
681        cnew.fMultiplicity=nstripsP;
682        cnew.fMultiplicityN=nstripsN;
683        cnew.fQErr=signalCut;
684        cnew.fNtracks=NTracks;
685
686        fITS->AddCluster(2,&cnew);
687
688        AliITSRecPoint rnew;
689        rnew.SetX(P*kconv);
690        rnew.SetZ(N*kconv);
691        rnew.SetQ(signal);
692        rnew.SetdEdX(dedx);
693        rnew.SetSigmaX2( kRMSx* kRMSx); 
694        rnew.SetSigmaZ2( kRMSz* kRMSz);
695        rnew.fTracks[0]=tr[0];
696        rnew.fTracks[1]=tr[1];
697        rnew.fTracks[2]=tr[2];
698
699        fITS->AddRecPoint(rnew);
700
701        return kTRUE;
702      }
703   }
704      return kFALSE;  
705 }
706
707
708 //------------------------------------------------------
709 void  AliITSClusterFinderSSD::CalcStepFactor(Float_t Psteo, Float_t Nsteo )
710 {
711 // calculate the step factor for matching clusters
712
713
714   // 95 is the pitch, 4000 - dimension along z ?
715
716
717   Float_t dz=fSegmentation->Dz();
718
719   fSFF = ( (Int_t)  (Psteo*dz/fPitch ) );// +1;
720   fSFB = ( (Int_t)  (Nsteo*dz/fPitch ) );// +1;
721
722 }
723
724
725 //-----------------------------------------------------------
726 AliITSclusterSSD* AliITSClusterFinderSSD::GetPSideCluster(Int_t idx)
727 {
728 // get P side clusters
729
730
731
732
733   if((idx<0)||(idx>=fNClusterP))
734     {
735       printf("AliITSClusterFinderSSD::GetPSideCluster  : index out of range\n");
736       return 0;
737     }
738   else
739     {
740       return (AliITSclusterSSD*)((*fClusterP)[idx]);
741     }
742 }
743
744 //-------------------------------------------------------
745 AliITSclusterSSD* AliITSClusterFinderSSD::GetNSideCluster(Int_t idx)
746 {
747 // get N side clusters
748
749   if((idx<0)||(idx>=fNClusterN))
750     {
751       printf("AliITSClusterFinderSSD::GetNSideCluster  : index out of range\n");
752       return 0;
753     }
754   else
755     {
756       return (AliITSclusterSSD*)((*fClusterN)[idx]);
757     }
758 }
759
760 //--------------------------------------------------------
761 AliITSclusterSSD* AliITSClusterFinderSSD::GetCluster(Int_t idx, Bool_t side)
762 {
763 // Get cluster
764
765   return (side) ? GetPSideCluster(idx) : GetNSideCluster(idx);
766 }
767
768 //_______________________________________________________________________
769
770 Bool_t AliITSClusterFinderSSD::GetCrossing (Float_t &P, Float_t &N) 
771
772 // get crossing 
773
774   Float_t Dx = fSegmentation->Dx(); // detector size in x direction, microns
775   Float_t Dz = fSegmentation->Dz(); // detector size in z direction, microns
776
777   Float_t xL; // x local coordinate
778   Float_t zL; // z local coordinate
779   Float_t x;  // x = xL + Dx/2
780   Float_t z;  // z = zL + Dz/2
781   Float_t xP; // x coordinate in the P side from the first P strip
782   Float_t xN; // x coordinate in the N side from the first N strip
783   Float_t StereoP,StereoN;
784
785   fSegmentation->Angles(StereoP,StereoN);
786   fTanP=TMath::Tan(StereoP);
787   fTanN=TMath::Tan(StereoN);
788   Float_t kP = fTanP; // Tangent of 0.0075 mrad
789   Float_t kN = fTanN; // Tangent of 0.0275 mrad
790
791     P *= fPitch;
792     N *= fPitch; 
793     xP = N;      // change the mistake for the P/N
794     xN = P;      // coordinates correspondence in this function
795     x = xP + kP*(Dz*kN-xP+xN)/(kP+kN);
796     z = (Dz*kN-xP+xN)/(kP+kN); 
797     xL = x - Dx/2;
798     zL = z - Dz/2;
799     P = xL;
800     N = zL;  
801
802     if(TMath::Abs(xL) > Dx/2 || TMath::Abs(zL) > Dz/2) return kFALSE;
803     // Check that xL and zL are inside the detector for the 
804     // correspondent xP and xN coordinates
805
806     return kTRUE;   
807 }
808
809
810 //_________________________________________________________________________
811
812 void AliITSClusterFinderSSD::GetCrossingError(Float_t& dP, Float_t& dN)
813 {
814 // get crossing error
815
816   Float_t dz, dx;
817   
818   dz = TMath::Abs(( dP + dN )*fPitch/(fTanP + fTanN) );
819   dx = fPitch*(TMath::Abs(dP*(1 - fTanP/(fTanP + fTanN))) +
820                TMath::Abs(dN *fTanP/(fTanP + fTanN) ));
821   
822   dN = dz;
823   dP = dx;
824 }
825
826
827
828
829
830
831