RC12, RC17 violation: suppression
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderV2SSD.cxx
1 /**************************************************************************
2  * Copyright(c) 2007-2009, 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 /* $Id$ */
17
18 ////////////////////////////////////////////////////////////////////////////
19 //            Implementation of the ITS clusterer V2 class                //
20 //                                                                        //
21 //          Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch            //
22 //          Last revision: 13-05-09 Enrico Fragiacomo                     //
23 //                                  enrico.fragiacomo@ts.infn.it          //
24 //                                                                        //
25 ///////////////////////////////////////////////////////////////////////////
26
27 #include <Riostream.h>
28 #include "AliLog.h"
29
30 #include "AliITSClusterFinderV2SSD.h"
31 #include "AliITSRecPoint.h"
32 #include "AliITSgeomTGeo.h"
33 #include "AliITSDetTypeRec.h"
34 #include "AliRawReader.h"
35 #include "AliITSRawStreamSSD.h"
36 #include <TClonesArray.h>
37 #include "AliITSdigitSSD.h"
38 #include "AliITSReconstructor.h"
39 #include "AliITSCalibrationSSD.h"
40 #include "AliITSsegmentationSSD.h"
41
42 Short_t *AliITSClusterFinderV2SSD::fgPairs = 0x0;
43 Int_t    AliITSClusterFinderV2SSD::fgPairsSize = 0;
44 const Float_t  AliITSClusterFinderV2SSD::fgkThreshold = 5.;
45
46 const Float_t AliITSClusterFinderV2SSD::fgkCosmic2008StripShifts[16][9] = 
47   {{-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35},  // DDL 512
48    {-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35,-0.35},  // DDL 513
49    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15},  // DDL 514
50    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15},  // DDL 515
51    { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00},  // DDL 516
52    { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00},  // DDL 517
53    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15},  // DDL 518
54    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15},  // DDL 519
55    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.25,-0.15},  // DDL 520
56    {-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15,-0.15},  // DDL 521
57    {-0.10,-0.10,-0.10,-0.40,-0.40,-0.40,-0.10,-0.10,-0.45},  // DDL 522
58    {-0.10,-0.10,-0.10,-0.35,-0.35,-0.35,-0.10,-0.35,-0.50},  // DDL 523
59    { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00},  // DDL 524
60    { 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00},  // DDL 525
61    { 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35, 0.35},  // DDL 526
62    { 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45, 0.45}}; // DDL 527
63
64 ClassImp(AliITSClusterFinderV2SSD)
65
66
67 AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),
68                                                                              fLastSSD1(AliITSgeomTGeo::GetModuleIndex(6,1,1)-1)
69 {
70 //Default constructor
71
72 }
73  
74 //______________________________________________________________________
75 AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(const AliITSClusterFinderV2SSD &cf) : AliITSClusterFinder(cf),                                               fLastSSD1(cf.fLastSSD1)
76 {
77   // Copy constructor
78 }
79
80 //______________________________________________________________________
81 AliITSClusterFinderV2SSD& AliITSClusterFinderV2SSD::operator=(const AliITSClusterFinderV2SSD&  cf ){
82   // Assignment operator
83
84   this->~AliITSClusterFinderV2SSD();
85   new(this) AliITSClusterFinderV2SSD(cf);
86   return *this;
87 }
88
89
90 void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){
91
92   //Find clusters V2
93   SetModule(mod);
94   FindClustersSSD(fDigits);
95
96 }
97
98 void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) {
99   //------------------------------------------------------------
100   // Actual SSD cluster finder
101   //------------------------------------------------------------
102   Int_t smaxall=alldigits->GetEntriesFast();
103   if (smaxall==0) return;
104
105
106   //---------------------------------------
107   // load recoparam and calibration
108   // 
109   static AliITSRecoParam *repa = NULL;  
110   if(!repa){
111     repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
112     if(!repa){
113       repa = AliITSRecoParam::GetHighFluxParam();
114       AliWarning("Using default AliITSRecoParam class");
115     }
116   }
117
118   AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
119   Float_t gain=0;
120   Float_t noise=0;
121   //---------------------------------------
122
123
124   //------------------------------------
125   // fill the digits array with zero-suppression condition
126   // Signal is converted in KeV
127   //
128   TObjArray digits;
129   for (Int_t i=0;i<smaxall; i++){
130     AliITSdigitSSD *d=(AliITSdigitSSD*)alldigits->UncheckedAt(i);
131
132     if(d->IsSideP()) noise = cal->GetNoiseP(d->GetStripNumber());  
133     else noise = cal->GetNoiseN(d->GetStripNumber());
134     if (d->GetSignal()<3.*noise) continue;
135
136     if(d->IsSideP()) gain = cal->GetGainP(d->GetStripNumber());  
137     else gain = cal->GetGainN(d->GetStripNumber());
138
139     Float_t q=gain*d->GetSignal(); //
140     q=cal->ADCToKeV(q); // converts the charge in KeV from ADC units
141     d->SetSignal(Int_t(q));
142
143     digits.AddLast(d);
144   }
145   Int_t smax = digits.GetEntriesFast();
146   if (smax==0) return;
147   //------------------------------------
148
149   
150   const Int_t kMax=1000;
151   Int_t np=0, nn=0; 
152   Ali1Dcluster pos[kMax], neg[kMax];
153   Float_t y=0., q=0., qmax=0.; 
154   Int_t lab[4]={-2,-2,-2,-2};  
155   Bool_t flag5 = 0;
156   
157   /*
158   cout<<"-----------------------------"<<endl;
159   cout<<"this is module "<<fModule;
160   cout<<endl;
161   cout<<endl;
162   */
163
164   //--------------------------------------------------------
165   // start 1D-clustering from the first digit in the digits array
166   //
167   AliITSdigitSSD *d=(AliITSdigitSSD*)digits.UncheckedAt(0);
168   q += d->GetSignal();
169   y += d->GetCoord2()*d->GetSignal();
170   qmax=d->GetSignal();
171   lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
172   
173   if(d->IsSideP()) {
174     noise = cal->GetNoiseP(d->GetStripNumber());  
175     gain = cal->GetGainP(d->GetStripNumber());
176   }
177   else {
178     noise = cal->GetNoiseN(d->GetStripNumber());
179     gain = cal->GetGainN(d->GetStripNumber());
180   }  
181   noise*=gain;
182   noise=cal->ADCToKeV(noise); // converts noise in KeV from ADC units
183
184   if(qmax>fgkThreshold*noise) flag5=1; // seed for the cluster
185
186   /*
187   cout<<d->GetSignal()<<" "<<noise<<" "<<flag5<<" "<<
188     d->GetCoord1()<<" "<<d->GetCoord2()<<endl;
189   */
190
191   Int_t curr=d->GetCoord2();
192   Int_t flag=d->GetCoord1();
193
194   // Note: the first side which will be processed is supposed to be the 
195   // P-side which is neg
196   Int_t *n=&nn;
197   Ali1Dcluster *c=neg;
198   if(flag) {n=&np; c=pos;} // in case we have only Nstrips (P was bad!)
199
200   Int_t nd=1;
201   Int_t milab[10];
202   for (Int_t ilab=0;ilab<10;ilab++){
203     milab[ilab]=-2;
204   }
205   milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2);
206
207
208   //----------------------------------------------------------
209   // search for neighboring digits
210   //
211   for (Int_t s=1; s<smax; s++) {
212       d=(AliITSdigitSSD*)digits.UncheckedAt(s);      
213       Int_t strip=d->GetCoord2();
214
215       // if digits is not a neighbour or side did not change 
216       // and at least one of the previous digits met the seed condition
217       // then creates a new 1D cluster
218       if ( ( ((strip-curr) > 1) || (flag!=d->GetCoord1()) ) ) {
219
220         if(flag5) {
221           //cout<<"here1"<<endl;
222          c[*n].SetY(y/q);
223          c[*n].SetQ(q);
224          c[*n].SetNd(nd);
225          CheckLabels2(milab);
226          c[*n].SetLabels(milab);
227
228          if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
229            // Note: fUseUnfoldingInClusterFinderSSD=kFALSE by default in RecoParam
230
231            //Split suspiciously big cluster
232            if (nd>4&&nd<25) {
233              c[*n].SetY(y/q-0.25*nd);
234              c[*n].SetQ(0.5*q);
235              (*n)++;
236              if (*n==kMax) {
237                Error("FindClustersSSD","Too many 1D clusters !");
238                return;
239              }
240              c[*n].SetY(y/q+0.25*nd);
241              c[*n].SetQ(0.5*q);
242              c[*n].SetNd(nd);
243              c[*n].SetLabels(milab);
244            }     
245            
246          } // unfolding is on
247
248          (*n)++;
249          if (*n==kMax) {
250           Error("FindClustersSSD","Too many 1D clusters !");
251           return;
252          }
253
254         } // flag5 set
255
256          // reset everything
257          y=q=qmax=0.;
258          nd=0;
259          flag5=0;
260          lab[0]=lab[1]=lab[2]=-2;
261          for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
262
263          // if side changed from P to N, switch to pos 1D clusters
264          // (if for some reason the side changed from N to P then do the opposite)
265          if (flag!=d->GetCoord1()) 
266            { if(!flag) {n=&np; c=pos;} else {n=&nn; c=neg;} }
267
268       } // end create new 1D cluster from previous neighboring digits
269
270       // continues adding digits to the previous cluster 
271       // or start a new one 
272       flag=d->GetCoord1();
273       q += d->GetSignal();
274       y += d->GetCoord2()*d->GetSignal();
275       nd++;
276
277       if(d->IsSideP()) {
278         noise = cal->GetNoiseP(d->GetStripNumber());  
279         gain = cal->GetGainP(d->GetStripNumber());
280       }
281       else {
282         noise = cal->GetNoiseN(d->GetStripNumber());
283         gain = cal->GetGainN(d->GetStripNumber());
284       }
285       noise*=gain;
286       noise=cal->ADCToKeV(noise); // converts the charge in KeV from ADC units
287
288       if(d->GetSignal()>fgkThreshold*noise) flag5=1;
289
290       /*
291   cout<<d->GetSignal()<<" "<<noise<<" "<<flag5<<" "<<
292     d->GetCoord1()<<" "<<d->GetCoord2()<<endl;
293       */
294
295       if (d->GetSignal()>qmax) {
296          qmax=d->GetSignal();
297          lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
298       }
299       for (Int_t ilab=0;ilab<10;ilab++) {
300         if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab))); 
301       }
302       curr=strip;
303
304
305   } // loop over digits, no more digits in the digits array
306
307
308   // add the last 1D cluster 
309   if(flag5) {
310
311     //  cout<<"here2"<<endl;
312
313     c[*n].SetY(y/q);
314     c[*n].SetQ(q);
315     c[*n].SetNd(nd);
316     c[*n].SetLabels(lab);
317     
318     if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
319       
320       //Split suspiciously big cluster
321       if (nd>4 && nd<25) {
322         c[*n].SetY(y/q-0.25*nd);
323         c[*n].SetQ(0.5*q);
324         (*n)++;
325         if (*n==kMax) {
326           Error("FindClustersSSD","Too many 1D clusters !");
327           return;
328         }
329         c[*n].SetY(y/q+0.25*nd);
330         c[*n].SetQ(0.5*q);
331         c[*n].SetNd(nd);
332         c[*n].SetLabels(lab);
333       }
334     } // unfolding is on
335     
336     (*n)++;
337     if (*n==kMax) {
338       Error("FindClustersSSD","Too many 1D clusters !");
339       return;
340     }
341
342   } // if flag5 last 1D cluster added
343
344
345   //------------------------------------------------------
346   // call FindClustersSSD to pair neg and pos 1D clusters 
347   // and create recpoints from the crosses
348   // Note1: neg are Pside and pos are Nside!!
349   // Note2: if there are no Pside digits nn=0 (bad strips!!) (same for Nside)
350   //
351   //  cout<<nn<<" Pside and "<<np<<" Nside clusters"<<endl;
352   FindClustersSSD(neg, nn, pos, np);
353   //
354   //-----------------------------------------------------
355
356 }
357
358
359 void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){
360
361     //------------------------------------------------------------
362   // This function creates ITS clusters from raw data
363   //------------------------------------------------------------
364   rawReader->Reset();
365   AliITSRawStreamSSD inputSSD(rawReader);
366   FindClustersSSD(&inputSSD,clusters);
367   
368 }
369
370 void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStreamSSD* input, 
371                                         TClonesArray** clusters) 
372 {
373   //------------------------------------------------------------
374   // Actual SSD cluster finder for raw data
375   //------------------------------------------------------------
376
377   static AliITSRecoParam *repa = NULL;
378   if(!repa){
379     repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
380     if(!repa){
381       repa = AliITSRecoParam::GetHighFluxParam();
382       AliWarning("Using default AliITSRecoParam class");
383     }
384   }
385
386   Int_t nClustersSSD = 0;
387   const Int_t kMax = 1000;
388   Ali1Dcluster clusters1D[2][kMax];
389   Int_t nClusters[2] = {0, 0};
390   Int_t lab[3]={-2,-2,-2};
391   Float_t q = 0.;
392   Float_t y = 0.;
393   Int_t nDigits = 0;
394   Float_t gain=0;
395   Float_t noise=0.;
396   //  Float_t pedestal=0.;
397   Float_t oldnoise=0.;
398   AliITSCalibrationSSD* cal=NULL;
399
400   Int_t matrix[12][1536];
401   Int_t iddl=-1;
402   Int_t iad=-1;
403   Int_t oddl = -1;
404   Int_t oad = -1;
405   Int_t oadc = -1;
406   Int_t ostrip = -1;
407   Int_t osignal = 65535;
408   Int_t n=0;
409   Bool_t next=0;
410
411   // read raw data input stream
412   while (kTRUE) {
413     
414     // reset signal matrix
415     for(Int_t i=0; i<12; i++) { for(Int_t j=0; j<1536; j++) { matrix[i][j] = 65535;} }
416     
417     if((osignal!=65535)&&(ostrip>0)&&(ostrip<1536)) { 
418       n++;
419       matrix[oadc][ostrip] = osignal; // recover data from previous occurence of input->Next() 
420     }
421     
422     // buffer data for ddl=iddl and ad=iad
423     while(kTRUE) {
424         
425       next = input->Next();
426       if((!next)&&(input->flag)) continue;
427       Int_t ddl=input->GetDDL(); 
428       Int_t ad=input->GetAD();
429       Int_t adc = input->GetADC(); adc = (adc<6)? adc : adc - 2;
430       Int_t strip = input->GetStrip();
431       if(input->GetSideFlag()) strip=1535-strip;
432       Int_t signal = input->GetSignal();
433
434       if((ddl==iddl)&&(ad==iad)&&(strip>0)&&(strip<1536)) {n++; matrix[adc][strip] = signal;}
435       else {if ((strip<1536) && (strip>0)) {oddl=iddl; oad=iad; oadc = adc; ostrip = strip; osignal=signal; iddl=ddl; iad=ad; break;}}
436       
437       if(!next)  {oddl=iddl; oad=iad; oadc = adc; ostrip = strip; osignal=signal; iddl=ddl; iad=ad; break;}
438       //break;
439     }
440     
441     // No SSD data
442     if(!next && oddl<0) break;
443     
444     if(n==0) continue; // first occurence
445     n=0; //osignal=0;
446
447     Float_t dStrip = 0;
448     if (repa->GetUseCosmicRunShiftsSSD()) {  // Special condition for 2007/2008 cosmic data
449       dStrip = fgkCosmic2008StripShifts[oddl][oad-1];
450     }
451     if (TMath::Abs(dStrip) > 1.5)
452       AliError(Form("Indexing error ? oddl = %d, dStrip %f\n",oddl,dStrip));
453     // fill 1Dclusters
454     for(Int_t iadc=0; iadc<12; iadc++) {  // loop over ADC index for ddl=oddl and ad=oad
455       
456       Int_t iimod = (oad - 1)  * 12 + iadc;
457       Int_t iModule = AliITSRawStreamSSD::GetModuleNumber(oddl,iimod);
458       if(iModule==-1) continue;
459       cal = (AliITSCalibrationSSD*)GetResp(iModule);
460
461       Bool_t first = 0;
462       Bool_t flag5 = 0;
463       
464       /*
465       for(Int_t istrip=0; istrip<768; istrip++) { // P-side
466         Int_t signal = matrix[iadc][istrip];
467         pedestal = cal->GetPedestalP(istrip);
468         matrix[iadc][istrip]=signal-(Int_t)pedestal;
469       } 
470       */
471
472       /*
473       Float_t cmode=0;
474       for(Int_t l=0; l<6; l++) {
475         cmode=0;
476         for(Int_t n=20; n<108; n++) cmode+=matrix[iadc][l*128+n];
477         cmode/=88.;
478         for(Int_t n=0; n<128; n++) matrix[iadc][l*128+n]-=(Int_t)cmode;
479         
480       }
481       */
482
483       Int_t istrip=0;
484       for(istrip=0; istrip<768; istrip++) { // P-side
485         
486         Int_t signal = TMath::Abs(matrix[iadc][istrip]);
487         
488         oldnoise = noise;
489         noise = cal->GetNoiseP(istrip); if(noise<1.) signal = 65535;
490         if(signal<3*noise) signal = 65535; // in case ZS was not done in hw do it now
491
492         //        if(cal->IsPChannelBad(istrip)) signal=0;
493
494         if (signal!=65535) {
495           gain = cal->GetGainP(istrip);
496           signal = (Int_t) ( signal * gain ); // signal is corrected for gain
497           if(signal>fgkThreshold*noise) flag5=1;
498           signal = (Int_t) cal->ADCToKeV( signal ); // signal is  converted in KeV 
499           
500           q += signal;    // add digit to current cluster
501           y += istrip * signal;   
502           nDigits++;
503           first=1;
504         }
505         
506         else if(first) {
507           
508           if ( (nDigits>0) && flag5 ) {
509             
510             Ali1Dcluster& cluster = clusters1D[0][nClusters[0]++];
511
512             if(q!=0) cluster.SetY(y/q + dStrip);
513             else cluster.SetY(istrip + dStrip -1);
514
515             cluster.SetQ(q);
516             cluster.SetNd(nDigits);
517             cluster.SetLabels(lab);
518             
519             if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
520               
521               //Split suspiciously big cluster
522               if (nDigits > 4&&nDigits < 25) {
523                 if(q!=0) cluster.SetY(y/q + dStrip - 0.25*nDigits);
524                 else cluster.SetY(istrip-1 + dStrip - 0.25*nDigits);
525                 cluster.SetQ(0.5*q);
526                 if (nClusters[0] == kMax) {
527                   Error("FindClustersSSD", "Too many 1D clusters !");
528                   return;
529                 }
530                 Ali1Dcluster& cluster2 = clusters1D[0][nClusters[0]++];
531                 if(q!=0) cluster2.SetY(y/q + dStrip + 0.25*nDigits);
532                 else cluster2.SetY(istrip-1 + dStrip + 0.25*nDigits);
533                 cluster2.SetQ(0.5*q);
534                 cluster2.SetNd(nDigits);
535                 cluster2.SetLabels(lab);
536               }
537             } // unfolding is on            
538           }
539           
540           y = q = 0.;
541           nDigits = 0;
542           first=0;
543           flag5=0;
544         }
545         
546       } // loop over strip on P-side
547       
548       // if last strip does have signal
549       if(first) {
550         
551          if ( (nDigits>0) && flag5 ) {
552           
553           Ali1Dcluster& cluster = clusters1D[0][nClusters[0]++];
554
555           if(q!=0) cluster.SetY(y/q + dStrip);
556           else cluster.SetY(istrip - 1 + dStrip);
557
558           cluster.SetQ(q);
559           cluster.SetNd(nDigits);
560           cluster.SetLabels(lab);
561           
562           if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
563             
564             //Split suspiciously big cluster
565             if (nDigits > 4&&nDigits < 25) {
566               if(q!=0) cluster.SetY(y/q + dStrip - 0.25*nDigits);
567               else cluster.SetY(istrip-1 + dStrip - 0.25*nDigits);
568               cluster.SetQ(0.5*q);
569               if (nClusters[0] == kMax) {
570                 Error("FindClustersSSD", "Too many 1D clusters !");
571                 return;
572               }
573               Ali1Dcluster& cluster2 = clusters1D[0][nClusters[0]++];
574               if(q!=0) cluster2.SetY(y/q + dStrip + 0.25*nDigits);
575               else cluster2.SetY(istrip-1 + dStrip + 0.25*nDigits);
576               cluster2.SetQ(0.5*q);
577               cluster2.SetNd(nDigits);
578               cluster2.SetLabels(lab);
579             }
580           } // unfolding is on    
581           
582         }
583         y = q = 0.;
584         nDigits = 0;
585         first=0;
586         flag5=0;
587       }
588       
589       /*
590       for(Int_t istrip=768; istrip<1536; istrip++) { // P-side
591         Int_t signal = matrix[iadc][istrip];
592         pedestal = cal->GetPedestalN(1535-istrip);
593         matrix[iadc][istrip]=signal-(Int_t)pedestal;
594       } 
595       */
596
597       /*
598       for(Int_t l=6; l<12; l++) {
599         Float_t cmode=0;
600         for(Int_t n=20; n<108; n++) cmode+=matrix[iadc][l*128+n];
601         cmode/=88.;
602         for(Int_t n=0; n<128; n++) matrix[iadc][l*128+n]-=(Int_t)cmode;
603       }
604       */
605
606       oldnoise = 0.;
607       noise = 0.;
608       Int_t strip=0;
609       for(Int_t iistrip=768; iistrip<1536; iistrip++) { // N-side
610         
611         Int_t signal = TMath::Abs(matrix[iadc][iistrip]);
612         strip = 1535-iistrip;
613
614         oldnoise = noise;
615         noise = cal->GetNoiseN(strip); if(noise<1.) signal=65535;
616
617         //        if(cal->IsNChannelBad(strip)) signal=0;
618
619         if(signal<3*noise) signal = 65535; // in case ZS was not done in hw do it now
620
621         if (signal!=65535) {
622           gain = cal->GetGainN(strip);
623           signal = (Int_t) ( signal * gain); // signal is corrected for gain
624           if(signal>fgkThreshold*noise) flag5=1;
625           signal = (Int_t) cal->ADCToKeV( signal ); // signal is  converted in KeV 
626           
627           // add digit to current cluster
628           q += signal;
629           y += strip * signal;
630           nDigits++;
631           first=1;
632         }
633
634         else if(first) {
635           
636           if ( (nDigits>0) && flag5 ) {
637             
638             Ali1Dcluster& cluster = clusters1D[1][nClusters[1]++];
639
640             if(q!=0) cluster.SetY(y/q - dStrip);
641             else cluster.SetY(strip+1 - dStrip);
642
643             cluster.SetQ(q);
644             cluster.SetNd(nDigits);
645             cluster.SetLabels(lab);
646             
647             if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
648
649               //Split suspiciously big cluster
650               if (nDigits > 4&&nDigits < 25) {
651                 cluster.SetY(y/q - dStrip - 0.25*nDigits);
652                 cluster.SetQ(0.5*q);
653                 if (nClusters[1] == kMax) {
654                   Error("FindClustersSSD", "Too many 1D clusters !");
655                   return;
656                 }
657                 Ali1Dcluster& cluster2 = clusters1D[1][nClusters[1]++];
658                 cluster2.SetY(y/q - dStrip + 0.25*nDigits);
659                 cluster2.SetQ(0.5*q);
660                 cluster2.SetNd(nDigits);
661                 cluster2.SetLabels(lab);
662               }       
663             } // unfolding is on
664           } 
665
666           y = q = 0.;
667           nDigits = 0;
668           first=0;        
669         flag5=0;
670         }
671         
672       } // loop over strips on N-side
673
674       if(first) {
675         
676           if ( (nDigits>0) && flag5 ) {
677           
678           Ali1Dcluster& cluster = clusters1D[1][nClusters[1]++];
679           
680           if(q!=0) cluster.SetY(y/q - dStrip);
681           else cluster.SetY(strip - dStrip + 1);
682
683           cluster.SetQ(q);
684           cluster.SetNd(nDigits);
685           cluster.SetLabels(lab);
686           
687           if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
688             
689             //Split suspiciously big cluster
690             if (nDigits > 4&&nDigits < 25) {
691               if(q!=0) cluster.SetY(y/q - dStrip - 0.25*nDigits);
692               else cluster.SetY(strip+1 - dStrip - 0.25*nDigits);
693               cluster.SetQ(0.5*q);
694               if (nClusters[1] == kMax) {
695                 Error("FindClustersSSD", "Too many 1D clusters !");
696                 return;
697               }
698               Ali1Dcluster& cluster2 = clusters1D[1][nClusters[1]++];
699               if(q!=0) cluster2.SetY(y/q - dStrip + 0.25*nDigits);
700               else cluster2.SetY(strip+1 - dStrip + 0.25*nDigits);
701               cluster2.SetQ(0.5*q);
702               cluster2.SetNd(nDigits);
703               cluster2.SetLabels(lab);
704             }
705           } // unfolding is on      
706         }
707
708         y = q = 0.;
709         nDigits = 0;
710         first=0;          
711         flag5=0;
712       }
713       
714       // create recpoints
715       if((nClusters[0])&&(nClusters[1])) {
716         
717         clusters[iModule] = new TClonesArray("AliITSRecPoint");
718         fModule = iModule;
719         FindClustersSSD(&clusters1D[0][0], nClusters[0], 
720                         &clusters1D[1][0], nClusters[1], clusters[iModule]);
721         Int_t nClustersn = clusters[iModule]->GetEntriesFast();
722         nClustersSSD += nClustersn;
723       }
724
725       nClusters[0] = nClusters[1] = 0;
726       y = q = 0.;
727       nDigits = 0;
728
729     } // loop over adc
730
731     if(!next) break;
732   }
733   
734   Info("FindClustersSSD", "found clusters in ITS SSD: %d", nClustersSSD);
735 }
736
737 void AliITSClusterFinderV2SSD::
738 FindClustersSSD(Ali1Dcluster* neg, Int_t nn, 
739                 Ali1Dcluster* pos, Int_t np,
740                 TClonesArray *clusters) {
741   //------------------------------------------------------------
742   // Actual SSD cluster finder
743   //------------------------------------------------------------
744
745   const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule);
746
747   //---------------------------------------
748   // load recoparam
749   // 
750   static AliITSRecoParam *repa = NULL;  
751   if(!repa){
752     repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
753     if(!repa){
754       repa = AliITSRecoParam::GetHighFluxParam();
755       AliWarning("Using default AliITSRecoParam class");
756     }
757   }
758
759   TClonesArray &cl=*clusters;
760   
761   AliITSsegmentationSSD *seg = dynamic_cast<AliITSsegmentationSSD*>(fDetTypeRec->GetSegmentationModel(2));
762   if (fModule>fLastSSD1) 
763     seg->SetLayer(6);
764   else 
765     seg->SetLayer(5);
766
767   Float_t hwSSD = seg->Dx()*1e-4/2;
768   Float_t hlSSD = seg->Dz()*1e-4/2;
769
770   Int_t idet=fNdet[fModule];
771   Int_t ncl=0;
772
773   //
774   Int_t *cnegative = new Int_t[np];  
775   Int_t *cused1 = new Int_t[np];
776   Int_t *negativepair = new Int_t[10*np];
777   Int_t *cpositive = new Int_t[nn];  
778   Int_t *cused2 = new Int_t[nn];  
779   Int_t *positivepair = new Int_t[10*nn];  
780   for (Int_t i=0;i<np;i++) {cnegative[i]=0; cused1[i]=0;}
781   for (Int_t i=0;i<nn;i++) {cpositive[i]=0; cused2[i]=0;}
782   for (Int_t i=0;i<10*np;i++) {negativepair[i]=0;}
783   for (Int_t i=0;i<10*nn;i++) {positivepair[i]=0;}
784
785   if ((np*nn) > fgPairsSize) {
786
787     if (fgPairs) delete [] fgPairs;
788     fgPairsSize = 4*np*nn;
789     fgPairs = new Short_t[fgPairsSize];
790   }
791   memset(fgPairs,0,sizeof(Short_t)*np*nn);
792
793   //
794   // find available pairs
795   //
796   for (Int_t i=0; i<np; i++) {
797     Float_t yp=pos[i].GetY(); 
798     if ( (pos[i].GetQ()>0) && (pos[i].GetQ()<3) ) continue;
799     for (Int_t j=0; j<nn; j++) {
800       if ( (neg[j].GetQ()>0) && (neg[j].GetQ()<3) ) continue;
801       Float_t yn=neg[j].GetY();
802
803       Float_t xt, zt;
804       seg->GetPadCxz(yn, yp, xt, zt);
805       //cout<<yn<<" "<<yp<<" "<<xt<<" "<<zt<<endl;
806       
807       if (TMath::Abs(xt)<hwSSD)
808       if (TMath::Abs(zt)<hlSSD) {
809         Int_t in = i*10+cnegative[i];
810         Int_t ip = j*10+cpositive[j];
811         if ((in < 10*np) && (ip < 10*nn)) {
812           negativepair[in] =j;  //index
813           positivepair[ip] =i;
814           cnegative[i]++;  //counters
815           cpositive[j]++;       
816           fgPairs[i*nn+j]=100;
817         }
818         else
819           AliError(Form("Index out of range: ip=%d, in=%d",ip,in));
820       }
821     }
822   }
823
824   /* //
825   // try to recover points out of but close to the module boundaries 
826   //
827   for (Int_t i=0; i<np; i++) {
828     Float_t yp=pos[i].GetY(); 
829     if ( (pos[i].GetQ()>0) && (pos[i].GetQ()<3) ) continue;
830     for (Int_t j=0; j<nn; j++) {
831       if ( (neg[j].GetQ()>0) && (neg[j].GetQ()<3) ) continue;
832       // if both 1Dclusters have an other cross continue
833       if (cpositive[j]&&cnegative[i]) continue;
834       Float_t yn=neg[j].GetY();
835       
836       Float_t xt, zt;
837       seg->GetPadCxz(yn, yp, xt, zt);
838       
839       if (TMath::Abs(xt)<hwSSD+0.1)
840       if (TMath::Abs(zt)<hlSSD+0.15) {
841         // tag 1Dcluster (eventually will produce low quality recpoint)
842         if (cnegative[i]==0) pos[i].SetNd(100);  // not available pair
843         if (cpositive[j]==0) neg[j].SetNd(100);  // not available pair
844         Int_t in = i*10+cnegative[i];
845         Int_t ip = j*10+cpositive[j];
846         if ((in < 10*np) && (ip < 10*nn)) {
847           negativepair[in] =j;  //index
848           positivepair[ip] =i;
849           cnegative[i]++;  //counters
850           cpositive[j]++;       
851           fgPairs[i*nn+j]=100;
852         }
853         else
854           AliError(Form("Index out of range: ip=%d, in=%d",ip,in));
855       }
856     }
857   }
858   */
859
860   //
861   Float_t lp[6];
862   Int_t milab[10];
863   Double_t ratio;
864   
865
866   if(repa->GetUseChargeMatchingInClusterFinderSSD()==kTRUE) {
867
868
869     //
870     // sign gold tracks
871     //
872     for (Int_t ip=0;ip<np;ip++){
873       Float_t xbest=1000,zbest=1000,qbest=0;
874       //
875       // select gold clusters
876       if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){ 
877         Float_t yp=pos[ip].GetY(); 
878         Int_t j = negativepair[10*ip];      
879
880         if( (pos[ip].GetQ()==0) && (neg[j].GetQ() ==0) ) { 
881           // both bad, hence continue;    
882           // mark both as used (to avoid recover at the end)
883           cused1[ip]++; 
884           cused2[j]++;
885           continue;
886         }
887
888         ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
889         //cout<<"ratio="<<ratio<<endl;
890
891         // charge matching (note that if posQ or negQ is 0 -> ratio=1 and the following condition is met
892         if (TMath::Abs(ratio)>0.2) continue; // note: 0.2=3xsigma_ratio calculated in cosmics tests
893
894         //
895         Float_t yn=neg[j].GetY();
896         
897         Float_t xt, zt;
898         seg->GetPadCxz(yn, yp, xt, zt);
899         
900         xbest=xt; zbest=zt; 
901
902         
903         qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
904         if( (pos[ip].GetQ()==0)||(neg[j].GetQ()==0)) qbest*=2; // in case of bad strips on one side keep all charge from the other one
905         
906         {
907           Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
908           mT2L->MasterToLocal(loc,trk);
909           lp[0]=trk[1];
910           lp[1]=trk[2];
911         }
912         lp[4]=qbest;        //Q
913         for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
914         for (Int_t ilab=0;ilab<3;ilab++){
915           milab[ilab] = pos[ip].GetLabel(ilab);
916           milab[ilab+3] = neg[j].GetLabel(ilab);
917         }
918         //
919         CheckLabels2(milab);
920         milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
921         Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
922
923         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
924         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
925         // out-of-diagonal element of covariance matrix
926         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
927         else if ( (info[0]>1) && (info[1]>1) ) { 
928           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
929           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
930           lp[5]=-6.48e-05;
931         }
932         else {
933           lp[2]=4.80e-06;      // 0.00219*0.00219
934           lp[3]=0.0093;        // 0.0964*0.0964;
935           if (info[0]==1) {
936             lp[5]=-0.00014;
937           }
938           else { 
939             lp[2]=2.79e-06;    // 0.0017*0.0017; 
940             lp[3]=0.00935;     // 0.967*0.967;
941             lp[5]=-4.32e-05;
942           }
943         }
944
945         AliITSRecPoint * cl2;
946         
947         if(clusters){  // Note clusters != 0 when method is called for rawdata
948           
949           
950           cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
951           
952           cl2->SetChargeRatio(ratio);           
953           cl2->SetType(1);
954           fgPairs[ip*nn+j]=1;
955
956           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
957             cl2->SetType(2);
958             fgPairs[ip*nn+j]=2;
959           }
960
961           if(pos[ip].GetQ()==0) cl2->SetType(3);
962           if(neg[j].GetQ()==0) cl2->SetType(4);
963
964           cused1[ip]++;
965           cused2[j]++;
966           
967         }
968         else{ // Note clusters == 0 when method is called for digits
969           
970           cl2 = new AliITSRecPoint(milab,lp,info);      
971           
972           cl2->SetChargeRatio(ratio);           
973           cl2->SetType(1);
974           fgPairs[ip*nn+j]=1;
975
976           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
977             cl2->SetType(2);
978             fgPairs[ip*nn+j]=2;
979           }
980
981           if(pos[ip].GetQ()==0) cl2->SetType(3);
982           if(neg[j].GetQ()==0) cl2->SetType(4);
983
984           cused1[ip]++;
985           cused2[j]++;
986
987           fDetTypeRec->AddRecPoint(*cl2);
988         }
989         ncl++;
990       }
991     }
992     
993     for (Int_t ip=0;ip<np;ip++){
994       Float_t xbest=1000,zbest=1000,qbest=0;
995       //
996       //
997       // select "silber" cluster
998       if ( cnegative[ip]==1 && cpositive[negativepair[10*ip]]==2){
999         Int_t in  = negativepair[10*ip];
1000         Int_t ip2 = positivepair[10*in];
1001         if (ip2==ip) ip2 =  positivepair[10*in+1];
1002         Float_t pcharge = pos[ip].GetQ()+pos[ip2].GetQ();
1003         
1004
1005
1006         ratio = (pcharge-neg[in].GetQ())/(pcharge+neg[in].GetQ());
1007         if ( (TMath::Abs(ratio)<0.2) && (pcharge!=0) ) {
1008           //if ( (TMath::Abs(pcharge-neg[in].GetQ())<30) && (pcharge!=0) ) { // 
1009           
1010           //
1011           // add first pair
1012           if ( (fgPairs[ip*nn+in]==100)&&(pos[ip].GetQ() ) ) {  //
1013             
1014             Float_t yp=pos[ip].GetY(); 
1015             Float_t yn=neg[in].GetY();
1016             
1017             Float_t xt, zt;
1018             seg->GetPadCxz(yn, yp, xt, zt);
1019             
1020             xbest=xt; zbest=zt; 
1021
1022             qbest =pos[ip].GetQ();
1023             Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1024             mT2L->MasterToLocal(loc,trk);
1025             lp[0]=trk[1];
1026             lp[1]=trk[2];
1027             
1028             lp[4]=qbest;        //Q
1029             for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1030             for (Int_t ilab=0;ilab<3;ilab++){
1031               milab[ilab] = pos[ip].GetLabel(ilab);
1032               milab[ilab+3] = neg[in].GetLabel(ilab);
1033             }
1034             //
1035             CheckLabels2(milab);
1036             ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ());
1037             milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det
1038             Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
1039             
1040         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1041         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1042         // out-of-diagonal element of covariance matrix
1043         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1044         else if ( (info[0]>1) && (info[1]>1) ) { 
1045           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1046           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1047           lp[5]=-6.48e-05;
1048         }
1049         else {
1050           lp[2]=4.80e-06;      // 0.00219*0.00219
1051           lp[3]=0.0093;        // 0.0964*0.0964;
1052           if (info[0]==1) {
1053             lp[5]=-0.00014;
1054           }
1055           else { 
1056             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1057             lp[3]=0.00935;     // 0.967*0.967;
1058             lp[5]=-4.32e-05;
1059           }
1060         }
1061
1062         AliITSRecPoint * cl2;
1063             if(clusters){
1064               
1065               cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1066               cl2->SetChargeRatio(ratio);       
1067               cl2->SetType(5);
1068               fgPairs[ip*nn+in] = 5;
1069               if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
1070                 cl2->SetType(6);
1071                 fgPairs[ip*nn+in] = 6;
1072               }     
1073             }
1074             else{
1075               cl2 = new AliITSRecPoint(milab,lp,info);
1076               cl2->SetChargeRatio(ratio);       
1077               cl2->SetType(5);
1078               fgPairs[ip*nn+in] = 5;
1079               if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
1080                 cl2->SetType(6);
1081                 fgPairs[ip*nn+in] = 6;
1082               }
1083               
1084               fDetTypeRec->AddRecPoint(*cl2);
1085             }
1086             ncl++;
1087           }
1088           
1089           
1090           //
1091           // add second pair
1092           
1093           //    if (!(cused1[ip2] || cused2[in])){  //
1094           if ( (fgPairs[ip2*nn+in]==100) && (pos[ip2].GetQ()) ) {
1095             
1096             Float_t yp=pos[ip2].GetY();
1097             Float_t yn=neg[in].GetY();
1098             
1099             Float_t xt, zt;
1100             seg->GetPadCxz(yn, yp, xt, zt);
1101             
1102             xbest=xt; zbest=zt; 
1103
1104             qbest =pos[ip2].GetQ();
1105             
1106             Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1107             mT2L->MasterToLocal(loc,trk);
1108             lp[0]=trk[1];
1109             lp[1]=trk[2];
1110             
1111             lp[4]=qbest;        //Q
1112             for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1113             for (Int_t ilab=0;ilab<3;ilab++){
1114               milab[ilab] = pos[ip2].GetLabel(ilab);
1115               milab[ilab+3] = neg[in].GetLabel(ilab);
1116             }
1117             //
1118             CheckLabels2(milab);
1119             ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ());
1120             milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det
1121             Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
1122             
1123         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1124         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1125         // out-of-diagonal element of covariance matrix
1126         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1127         else if ( (info[0]>1) && (info[1]>1) ) { 
1128           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1129           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1130           lp[5]=-6.48e-05;
1131         }
1132         else {
1133           lp[2]=4.80e-06;      // 0.00219*0.00219
1134           lp[3]=0.0093;        // 0.0964*0.0964;
1135           if (info[0]==1) {
1136             lp[5]=-0.00014;
1137           }
1138           else { 
1139             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1140             lp[3]=0.00935;     // 0.967*0.967;
1141             lp[5]=-4.32e-05;
1142           }
1143         }
1144
1145             AliITSRecPoint * cl2;
1146             if(clusters){
1147               cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1148               
1149               cl2->SetChargeRatio(ratio);       
1150               cl2->SetType(5);
1151               fgPairs[ip2*nn+in] =5;
1152               if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
1153                 cl2->SetType(6);
1154                 fgPairs[ip2*nn+in] =6;
1155               }
1156             }
1157             else{
1158               cl2 = new AliITSRecPoint(milab,lp,info);
1159               cl2->SetChargeRatio(ratio);       
1160               cl2->SetType(5);
1161               fgPairs[ip2*nn+in] =5;
1162               if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
1163                 cl2->SetType(6);
1164                 fgPairs[ip2*nn+in] =6;
1165               }
1166               
1167               fDetTypeRec->AddRecPoint(*cl2);
1168             }
1169             ncl++;
1170           }
1171           
1172           cused1[ip]++;
1173           cused1[ip2]++;
1174           cused2[in]++;
1175
1176         } // charge matching condition
1177
1178       } // 2 Pside cross 1 Nside
1179     } // loop over Pside clusters
1180     
1181     
1182       
1183       //  
1184     for (Int_t jn=0;jn<nn;jn++){
1185       if (cused2[jn]) continue;
1186       Float_t xbest=1000,zbest=1000,qbest=0;
1187       // select "silber" cluster
1188       if ( cpositive[jn]==1 && cnegative[positivepair[10*jn]]==2){
1189         Int_t ip  = positivepair[10*jn];
1190         Int_t jn2 = negativepair[10*ip];
1191         if (jn2==jn) jn2 =  negativepair[10*ip+1];
1192         Float_t pcharge = neg[jn].GetQ()+neg[jn2].GetQ();
1193         //
1194         
1195
1196         ratio = (pcharge-pos[ip].GetQ())/(pcharge+pos[ip].GetQ());
1197         if ( (TMath::Abs(ratio)<0.2) && (pcharge!=0) ) {
1198
1199           /*
1200         if ( (TMath::Abs(pcharge-pos[ip].GetQ())<30) &&  // charge matching 
1201              (pcharge!=0) ) { // reject combinations of bad strips
1202           */
1203
1204
1205           //
1206           // add first pair
1207           //    if (!(cused1[ip]||cused2[jn])){
1208           if ( (fgPairs[ip*nn+jn]==100) && (neg[jn].GetQ()) ) {  //
1209             
1210             Float_t yn=neg[jn].GetY(); 
1211             Float_t yp=pos[ip].GetY();
1212
1213             Float_t xt, zt;
1214             seg->GetPadCxz(yn, yp, xt, zt);
1215             
1216             xbest=xt; zbest=zt; 
1217
1218             qbest =neg[jn].GetQ();
1219
1220             {
1221               Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1222               mT2L->MasterToLocal(loc,trk);
1223               lp[0]=trk[1];
1224               lp[1]=trk[2];
1225           }
1226           
1227           lp[4]=qbest;        //Q
1228           for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1229           for (Int_t ilab=0;ilab<3;ilab++){
1230             milab[ilab] = pos[ip].GetLabel(ilab);
1231             milab[ilab+3] = neg[jn].GetLabel(ilab);
1232           }
1233           //
1234           CheckLabels2(milab);
1235           ratio = (pos[ip].GetQ()-neg[jn].GetQ())/(pos[ip].GetQ()+neg[jn].GetQ());
1236           milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det
1237           Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]};
1238
1239         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1240         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1241         // out-of-diagonal element of covariance matrix
1242         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1243         else if ( (info[0]>1) && (info[1]>1) ) { 
1244           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1245           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1246           lp[5]=-6.48e-05;
1247         }
1248         else {
1249           lp[2]=4.80e-06;      // 0.00219*0.00219
1250           lp[3]=0.0093;        // 0.0964*0.0964;
1251           if (info[0]==1) {
1252             lp[5]=-0.00014;
1253           }
1254           else { 
1255             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1256             lp[3]=0.00935;     // 0.967*0.967;
1257             lp[5]=-4.32e-05;
1258           }
1259         }
1260
1261           AliITSRecPoint * cl2;
1262           if(clusters){
1263             cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1264
1265             cl2->SetChargeRatio(ratio);         
1266             cl2->SetType(7);
1267             fgPairs[ip*nn+jn] =7;
1268             if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
1269               cl2->SetType(8);
1270               fgPairs[ip*nn+jn]=8;
1271             }
1272
1273           }
1274           else{
1275             cl2 = new AliITSRecPoint(milab,lp,info);
1276             cl2->SetChargeRatio(ratio);         
1277             cl2->SetType(7);
1278             fgPairs[ip*nn+jn] =7;
1279             if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
1280               cl2->SetType(8);
1281               fgPairs[ip*nn+jn]=8;
1282             }
1283
1284             fDetTypeRec->AddRecPoint(*cl2);
1285           }
1286           ncl++;
1287         }
1288         //
1289         // add second pair
1290         //      if (!(cused1[ip]||cused2[jn2])){
1291         if ( (fgPairs[ip*nn+jn2]==100)&&(neg[jn2].GetQ() ) ) {  //
1292
1293           Float_t yn=neg[jn2].GetY(); 
1294           Double_t yp=pos[ip].GetY(); 
1295
1296           Float_t xt, zt;
1297           seg->GetPadCxz(yn, yp, xt, zt);
1298           
1299           xbest=xt; zbest=zt; 
1300
1301           qbest =neg[jn2].GetQ();
1302
1303           {
1304           Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1305           mT2L->MasterToLocal(loc,trk);
1306           lp[0]=trk[1];
1307           lp[1]=trk[2];
1308           }
1309
1310           lp[4]=qbest;        //Q
1311           for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1312           for (Int_t ilab=0;ilab<3;ilab++){
1313             milab[ilab] = pos[ip].GetLabel(ilab);
1314             milab[ilab+3] = neg[jn2].GetLabel(ilab);
1315           }
1316           //
1317           CheckLabels2(milab);
1318           ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ());
1319           milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det
1320           Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]};
1321
1322         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1323         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1324         // out-of-diagonal element of covariance matrix
1325         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1326         else if ( (info[0]>1) && (info[1]>1) ) { 
1327           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1328           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1329           lp[5]=-6.48e-05;
1330         }
1331         else {
1332           lp[2]=4.80e-06;      // 0.00219*0.00219
1333           lp[3]=0.0093;        // 0.0964*0.0964;
1334           if (info[0]==1) {
1335             lp[5]=-0.00014;
1336           }
1337           else { 
1338             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1339             lp[3]=0.00935;     // 0.967*0.967;
1340             lp[5]=-4.32e-05;
1341           }
1342         }
1343
1344           AliITSRecPoint * cl2;
1345           if(clusters){
1346             cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1347
1348
1349             cl2->SetChargeRatio(ratio);         
1350             fgPairs[ip*nn+jn2]=7;
1351             cl2->SetType(7);
1352             if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
1353               cl2->SetType(8);
1354               fgPairs[ip*nn+jn2]=8;
1355             }
1356             
1357           }
1358           else{
1359             cl2 = new AliITSRecPoint(milab,lp,info);
1360             cl2->SetChargeRatio(ratio);         
1361             fgPairs[ip*nn+jn2]=7;
1362             cl2->SetType(7);
1363             if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
1364               cl2->SetType(8);
1365               fgPairs[ip*nn+jn2]=8;
1366             }
1367
1368             fDetTypeRec->AddRecPoint(*cl2);
1369           }
1370
1371           ncl++;
1372         }
1373         cused1[ip]++;
1374         cused2[jn]++;
1375         cused2[jn2]++;
1376
1377         } // charge matching condition
1378
1379       } // 2 Nside cross 1 Pside
1380     } // loop over Pside clusters
1381
1382   
1383     
1384     for (Int_t ip=0;ip<np;ip++){
1385
1386       if(cused1[ip]) continue;
1387
1388
1389       Float_t xbest=1000,zbest=1000,qbest=0;
1390       //
1391       // 2x2 clusters
1392       //
1393       if ( (cnegative[ip]==2) && cpositive[negativepair[10*ip]]==2){ 
1394         Float_t minchargediff =4.;
1395         Float_t minchargeratio =0.2;
1396
1397         Int_t j=-1;
1398         for (Int_t di=0;di<cnegative[ip];di++){
1399           Int_t   jc = negativepair[ip*10+di];
1400           Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1401           ratio = (pos[ip].GetQ()-neg[jc].GetQ())/(pos[ip].GetQ()+neg[jc].GetQ()); 
1402           //if (TMath::Abs(chargedif)<minchargediff){
1403           if (TMath::Abs(ratio)<0.2){
1404             j =jc;
1405             minchargediff = TMath::Abs(chargedif);
1406             minchargeratio = TMath::Abs(ratio);
1407           }
1408         }
1409         if (j<0) continue;  // not proper cluster      
1410         
1411
1412         Int_t count =0;
1413         for (Int_t di=0;di<cnegative[ip];di++){
1414           Int_t   jc = negativepair[ip*10+di];
1415           Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1416           if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1417         }
1418         if (count>1) continue;  // more than one "proper" cluster for positive
1419         //
1420         
1421         count =0;
1422         for (Int_t dj=0;dj<cpositive[j];dj++){
1423           Int_t   ic  = positivepair[j*10+dj];
1424           Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
1425           if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1426         }
1427         if (count>1) continue;  // more than one "proper" cluster for negative
1428         
1429         Int_t jp = 0;
1430         
1431         count =0;
1432         for (Int_t dj=0;dj<cnegative[jp];dj++){
1433           Int_t   ic = positivepair[jp*10+dj];
1434           Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
1435           if (TMath::Abs(chargedif)<minchargediff+4.) count++;
1436         }
1437         if (count>1) continue;   
1438         if (fgPairs[ip*nn+j]<100) continue;
1439         //
1440         
1441
1442
1443         //almost gold clusters
1444         Float_t yp=pos[ip].GetY(); 
1445         Float_t yn=neg[j].GetY();      
1446         Float_t xt, zt;
1447         seg->GetPadCxz(yn, yp, xt, zt); 
1448         xbest=xt; zbest=zt; 
1449         qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
1450         {
1451           Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1452           mT2L->MasterToLocal(loc,trk);
1453           lp[0]=trk[1];
1454           lp[1]=trk[2];
1455         }
1456         lp[4]=qbest;        //Q
1457         for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1458         for (Int_t ilab=0;ilab<3;ilab++){
1459           milab[ilab] = pos[ip].GetLabel(ilab);
1460           milab[ilab+3] = neg[j].GetLabel(ilab);
1461         }
1462         //
1463         CheckLabels2(milab);
1464         if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!!
1465         ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
1466         milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
1467         Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
1468
1469         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1470         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1471         // out-of-diagonal element of covariance matrix
1472         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1473         else if ( (info[0]>1) && (info[1]>1) ) { 
1474           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1475           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1476           lp[5]=-6.48e-05;
1477         }
1478         else {
1479           lp[2]=4.80e-06;      // 0.00219*0.00219
1480           lp[3]=0.0093;        // 0.0964*0.0964;
1481           if (info[0]==1) {
1482             lp[5]=-0.00014;
1483           }
1484           else { 
1485             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1486             lp[3]=0.00935;     // 0.967*0.967;
1487             lp[5]=-4.32e-05;
1488           }
1489         }
1490
1491         AliITSRecPoint * cl2;
1492         if(clusters){
1493           cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1494                   
1495           cl2->SetChargeRatio(ratio);           
1496           cl2->SetType(10);
1497           fgPairs[ip*nn+j]=10;
1498           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1499             cl2->SetType(11);
1500             fgPairs[ip*nn+j]=11;
1501           }
1502           cused1[ip]++;
1503           cused2[j]++;      
1504         }
1505         else{
1506           cl2 = new AliITSRecPoint(milab,lp,info);
1507           cl2->SetChargeRatio(ratio);           
1508           cl2->SetType(10);
1509           fgPairs[ip*nn+j]=10;
1510           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1511             cl2->SetType(11);
1512             fgPairs[ip*nn+j]=11;
1513           }
1514           cused1[ip]++;
1515           cused2[j]++;      
1516           
1517           fDetTypeRec->AddRecPoint(*cl2);
1518         }      
1519         ncl++;
1520         
1521       } // 2X2
1522     } // loop over Pside 1Dclusters
1523
1524
1525
1526     for (Int_t ip=0;ip<np;ip++){
1527
1528       if(cused1[ip]) continue;
1529
1530
1531       Float_t xbest=1000,zbest=1000,qbest=0;
1532       //
1533       // manyxmany clusters
1534       //
1535       if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){ 
1536         Float_t minchargediff =4.;
1537         Int_t j=-1;
1538         for (Int_t di=0;di<cnegative[ip];di++){
1539           Int_t   jc = negativepair[ip*10+di];
1540           Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1541           if (TMath::Abs(chargedif)<minchargediff){
1542             j =jc;
1543             minchargediff = TMath::Abs(chargedif);
1544           }
1545         }
1546         if (j<0) continue;  // not proper cluster      
1547         
1548         Int_t count =0;
1549         for (Int_t di=0;di<cnegative[ip];di++){
1550           Int_t   jc = negativepair[ip*10+di];
1551           Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1552           if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1553         }
1554         if (count>1) continue;  // more than one "proper" cluster for positive
1555         //
1556         
1557         count =0;
1558         for (Int_t dj=0;dj<cpositive[j];dj++){
1559           Int_t   ic  = positivepair[j*10+dj];
1560           Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
1561           if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1562         }
1563         if (count>1) continue;  // more than one "proper" cluster for negative
1564         
1565         Int_t jp = 0;
1566         
1567         count =0;
1568         for (Int_t dj=0;dj<cnegative[jp];dj++){
1569           Int_t   ic = positivepair[jp*10+dj];
1570           Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
1571           if (TMath::Abs(chargedif)<minchargediff+4.) count++;
1572         }
1573         if (count>1) continue;   
1574         if (fgPairs[ip*nn+j]<100) continue;
1575         //
1576         
1577         //almost gold clusters
1578         Float_t yp=pos[ip].GetY(); 
1579         Float_t yn=neg[j].GetY();
1580       
1581
1582         Float_t xt, zt;
1583         seg->GetPadCxz(yn, yp, xt, zt);
1584         
1585         xbest=xt; zbest=zt; 
1586
1587         qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
1588
1589         {
1590           Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1591           mT2L->MasterToLocal(loc,trk);
1592           lp[0]=trk[1];
1593           lp[1]=trk[2];
1594         }
1595         lp[4]=qbest;        //Q
1596         for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1597         for (Int_t ilab=0;ilab<3;ilab++){
1598           milab[ilab] = pos[ip].GetLabel(ilab);
1599           milab[ilab+3] = neg[j].GetLabel(ilab);
1600         }
1601         //
1602         CheckLabels2(milab);
1603         if ((neg[j].GetQ()==0)&&(pos[ip].GetQ()==0)) continue; // reject crosses of bad strips!!
1604         ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
1605         milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
1606         Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
1607
1608         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1609         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1610         // out-of-diagonal element of covariance matrix
1611         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1612         else if ( (info[0]>1) && (info[1]>1) ) { 
1613           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1614           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1615           lp[5]=-6.48e-05;
1616         }
1617         else {
1618           lp[2]=4.80e-06;      // 0.00219*0.00219
1619           lp[3]=0.0093;        // 0.0964*0.0964;
1620           if (info[0]==1) {
1621             lp[5]=-0.00014;
1622           }
1623           else { 
1624             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1625             lp[3]=0.00935;     // 0.967*0.967;
1626             lp[5]=-4.32e-05;
1627           }
1628         }
1629
1630         AliITSRecPoint * cl2;
1631         if(clusters){
1632           cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1633                   
1634           cl2->SetChargeRatio(ratio);           
1635           cl2->SetType(12);
1636           fgPairs[ip*nn+j]=12;
1637           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1638             cl2->SetType(13);
1639             fgPairs[ip*nn+j]=13;
1640           }
1641           cused1[ip]++;
1642           cused2[j]++;      
1643         }
1644         else{
1645           cl2 = new AliITSRecPoint(milab,lp,info);
1646           cl2->SetChargeRatio(ratio);           
1647           cl2->SetType(12);
1648           fgPairs[ip*nn+j]=12;
1649           if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1650             cl2->SetType(13);
1651             fgPairs[ip*nn+j]=13;
1652           }
1653           cused1[ip]++;
1654           cused2[j]++;      
1655           
1656           fDetTypeRec->AddRecPoint(*cl2);
1657         }      
1658         ncl++;
1659         
1660       } // manyXmany
1661     } // loop over Pside 1Dclusters
1662     
1663   } // use charge matching
1664   
1665   
1666   // recover all the other crosses
1667   //  
1668   for (Int_t i=0; i<np; i++) {
1669     Float_t xbest=1000,zbest=1000,qbest=0;
1670     Float_t yp=pos[i].GetY(); 
1671     if ((pos[i].GetQ()>0)&&(pos[i].GetQ()<3)) continue;
1672     for (Int_t j=0; j<nn; j++) {
1673     //    for (Int_t di = 0;di<cpositive[i];di++){
1674     //  Int_t j = negativepair[10*i+di];
1675       if ((neg[j].GetQ()>0)&&(neg[j].GetQ()<3)) continue;
1676
1677       if ((neg[j].GetQ()==0)&&(pos[i].GetQ()==0)) continue; // reject crosses of bad strips!!
1678
1679       if (cused2[j]||cused1[i]) continue;      
1680       if (fgPairs[i*nn+j]>0 &&fgPairs[i*nn+j]<100) continue;
1681       ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ());      
1682       Float_t yn=neg[j].GetY();
1683       
1684       Float_t xt, zt;
1685       seg->GetPadCxz(yn, yp, xt, zt);
1686       
1687       if (TMath::Abs(xt)<hwSSD)
1688       if (TMath::Abs(zt)<hlSSD) {
1689         xbest=xt; zbest=zt; 
1690
1691         qbest=0.5*(pos[i].GetQ()+neg[j].GetQ());
1692
1693         {
1694         Double_t loc[3]={xbest,0.,zbest},trk[3]={0.,0.,0.};
1695         mT2L->MasterToLocal(loc,trk);
1696         lp[0]=trk[1];
1697         lp[1]=trk[2];
1698         }
1699         lp[4]=qbest;        //Q
1700         for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1701         for (Int_t ilab=0;ilab<3;ilab++){
1702           milab[ilab] = pos[i].GetLabel(ilab);
1703           milab[ilab+3] = neg[j].GetLabel(ilab);
1704         }
1705         //
1706         CheckLabels2(milab);
1707         milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det
1708         Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
1709
1710         lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1711         lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1712         // out-of-diagonal element of covariance matrix
1713         if( (info[0]==1) && (info[1]==1) ) lp[5]=-0.00012;
1714         else if ( (info[0]>1) && (info[1]>1) ) { 
1715           lp[2]=2.63e-06;    // 0.0016*0.0016;  //SigmaY2
1716           lp[3]=0.0065;      // 0.08*0.08;   //SigmaZ2
1717           lp[5]=-6.48e-05;
1718         }
1719         else {
1720           lp[2]=4.80e-06;      // 0.00219*0.00219
1721           lp[3]=0.0093;        // 0.0964*0.0964;
1722           if (info[0]==1) {
1723             lp[5]=-0.00014;
1724           }
1725           else { 
1726             lp[2]=2.79e-06;    // 0.0017*0.0017; 
1727             lp[3]=0.00935;     // 0.967*0.967;
1728             lp[5]=-4.32e-05;
1729           }
1730         }
1731
1732         AliITSRecPoint * cl2;
1733         if(clusters){
1734           cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
1735
1736           cl2->SetChargeRatio(ratio);
1737           cl2->SetType(100+cpositive[j]+cnegative[i]);    
1738
1739           if(pos[i].GetQ()==0) cl2->SetType(200+cpositive[j]+cnegative[i]);
1740           if(neg[j].GetQ()==0) cl2->SetType(300+cpositive[j]+cnegative[i]);
1741
1742         }
1743         else{
1744           cl2 = new AliITSRecPoint(milab,lp,info);
1745           cl2->SetChargeRatio(ratio);
1746           cl2->SetType(100+cpositive[j]+cnegative[i]);
1747           
1748           if(pos[i].GetQ()==0) cl2->SetType(200+cpositive[j]+cnegative[i]);
1749           if(neg[j].GetQ()==0) cl2->SetType(300+cpositive[j]+cnegative[i]);
1750
1751           fDetTypeRec->AddRecPoint(*cl2);
1752         }
1753         ncl++;
1754       }
1755     }
1756   }
1757
1758
1759
1760   if(repa->GetUseBadChannelsInClusterFinderSSD()==kTRUE) {
1761     
1762     //---------------------------------------------------------
1763     // recover crosses of good 1D clusters with bad strips on the other side
1764     // Note1: at first iteration skip modules with a bad side (or almost), (would produce too many fake!) 
1765     // Note2: for modules with a bad side see below 
1766     
1767     AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
1768     Int_t countPbad=0, countNbad=0;
1769     for(Int_t ib=0; ib<768; ib++) {
1770       if(cal->IsPChannelBad(ib)) countPbad++;
1771       if(cal->IsNChannelBad(ib)) countNbad++;
1772     }
1773     //  AliInfo(Form("module %d has %d P- and %d N-bad strips",fModule,countPbad,countNbad));
1774     
1775     if( (countPbad<100) && (countNbad<100) ) { // no bad side!!
1776       
1777       for (Int_t i=0; i<np; i++) { // loop over Nside 1Dclusters with no crosses
1778         if(cnegative[i]) continue; // if intersecting Pside clusters continue;
1779         
1780         //      for(Int_t ib=0; ib<768; ib++) { // loop over all Pstrips
1781         for(Int_t ib=15; ib<753; ib++) { // loop over all Pstrips
1782           
1783           if(cal->IsPChannelBad(ib)) { // check if strips is bad
1784             Float_t yN=pos[i].GetY();   
1785             Float_t xt, zt;
1786             seg->GetPadCxz(1.*ib, yN, xt, zt);  
1787             
1788             //----------
1789             // bad Pstrip is crossing the Nside 1Dcluster -> create recpoint
1790             // 
1791             if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
1792               Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1793               mT2L->MasterToLocal(loc,trk);
1794               lp[0]=trk[1];
1795               lp[1]=trk[2];        
1796               lp[4]=pos[i].GetQ(); //Q
1797               for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1798               for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab);       
1799               CheckLabels2(milab);
1800               milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det
1801               Int_t info[3] = {pos[i].GetNd(),0,fNlayer[fModule]};
1802               
1803               lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1804               lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1805               lp[5]=-0.00012; // out-of-diagonal element of covariance matrix
1806               if (info[0]>1) {
1807                 lp[2]=4.80e-06;
1808                 lp[3]=0.0093;
1809                 lp[5]=0.00014;
1810               }
1811                       
1812               AliITSRecPoint * cl2;
1813               if(clusters){
1814                 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);          
1815                 cl2->SetChargeRatio(1.);
1816               cl2->SetType(50);   
1817               }
1818               else{
1819                 cl2 = new AliITSRecPoint(milab,lp,info);
1820                 cl2->SetChargeRatio(1.);
1821                 cl2->SetType(50);
1822                 fDetTypeRec->AddRecPoint(*cl2);
1823               }
1824               ncl++;
1825             } // cross is within the detector
1826             //
1827             //--------------
1828             
1829           } // bad Pstrip
1830           
1831         } // end loop over Pstrips
1832         
1833       } // end loop over Nside 1D clusters
1834       
1835       for (Int_t j=0; j<nn; j++) { // loop over Pside 1D clusters with no crosses
1836         if(cpositive[j]) continue;
1837         
1838         //      for(Int_t ib=0; ib<768; ib++) { // loop over all Nside strips
1839         for(Int_t ib=15; ib<753; ib++) { // loop over all Nside strips
1840           
1841           if(cal->IsNChannelBad(ib)) { // check if strip is bad
1842             Float_t yP=neg[j].GetY();   
1843             Float_t xt, zt;
1844             seg->GetPadCxz(yP, 1.*ib, xt, zt);  
1845             
1846             //----------
1847             // bad Nstrip is crossing the Pside 1Dcluster -> create recpoint
1848             // 
1849             if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
1850               Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1851               mT2L->MasterToLocal(loc,trk);
1852               lp[0]=trk[1];
1853               lp[1]=trk[2];        
1854               lp[4]=neg[j].GetQ(); //Q
1855               for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1856               for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab);       
1857               CheckLabels2(milab);
1858               milab[3]=( j << 10 ) + idet; // pos|neg|det
1859               Int_t info[3]={0,(Int_t)neg[j].GetNd(),fNlayer[fModule]};
1860               
1861               lp[2]=4.968e-06;     // 0.00223*0.00223;  //SigmaY2
1862               lp[3]=0.012;         // 0.110*0.110;  //SigmaZ2
1863               lp[5]=-0.00012; // out-of-diagonal element of covariance matrix
1864               if (info[0]>1) {
1865                 lp[2]=2.79e-06;
1866                 lp[3]=0.00935;
1867                 lp[5]=-4.32e-05;
1868               }
1869               
1870               AliITSRecPoint * cl2;
1871               if(clusters){
1872                 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);          
1873                 cl2->SetChargeRatio(1.);
1874                 cl2->SetType(60);         
1875               }
1876               else{
1877                 cl2 = new AliITSRecPoint(milab,lp,info);
1878                 cl2->SetChargeRatio(1.);
1879                 cl2->SetType(60);
1880                 fDetTypeRec->AddRecPoint(*cl2);
1881               }
1882               ncl++;
1883             } // cross is within the detector
1884             //
1885             //--------------
1886             
1887           } // bad Nstrip
1888         } // end loop over Nstrips
1889       } // end loop over Pside 1D clusters
1890       
1891     } // no bad sides 
1892     
1893     //---------------------------------------------------------
1894     
1895     else if( (countPbad>700) && (countNbad<100) ) { // bad Pside!!
1896       
1897       for (Int_t i=0; i<np; i++) { // loop over Nside 1Dclusters with no crosses
1898         if(cnegative[i]) continue; // if intersecting Pside clusters continue;
1899         
1900         Float_t xt, zt;
1901         Float_t yN=pos[i].GetY();       
1902         Float_t yP=0.;
1903         if (seg->GetLayer()==5) yP = yN + (7.6/1.9);
1904         else yP = yN - (7.6/1.9);
1905         seg->GetPadCxz(yP, yN, xt, zt); 
1906         
1907         if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
1908           Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1909           mT2L->MasterToLocal(loc,trk);
1910           lp[0]=trk[1];
1911           lp[1]=trk[2];        
1912           lp[4]=pos[i].GetQ(); //Q
1913           for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1914           for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = pos[i].GetLabel(ilab);   
1915           CheckLabels2(milab);
1916           milab[3]=( (i<<10) << 10 ) + idet; // pos|neg|det
1917           Int_t info[3] = {(Int_t)pos[i].GetNd(),0,fNlayer[fModule]};
1918           
1919           lp[2]=0.00098;    // 0.031*0.031;  //SigmaY2
1920           lp[3]=1.329;      // 1.15*1.15;  //SigmaZ2
1921           lp[5]=-0.0359;
1922           if(info[0]>1) lp[2]=0.00097;
1923
1924           AliITSRecPoint * cl2;
1925           if(clusters){
1926             cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);      
1927             cl2->SetChargeRatio(1.);
1928             cl2->SetType(70);     
1929           }
1930           else{
1931             cl2 = new AliITSRecPoint(milab,lp,info);
1932             cl2->SetChargeRatio(1.);
1933             cl2->SetType(70);
1934             fDetTypeRec->AddRecPoint(*cl2);
1935           }
1936           ncl++;
1937         } // cross is within the detector
1938         //
1939         //--------------
1940         
1941       } // end loop over Nside 1D clusters
1942       
1943     } // bad Pside module
1944     
1945     else if( (countNbad>700) && (countPbad<100) ) { // bad Nside!!
1946       
1947       for (Int_t j=0; j<nn; j++) { // loop over Pside 1D clusters with no crosses
1948         if(cpositive[j]) continue;
1949         
1950         Float_t xt, zt;
1951         Float_t yP=neg[j].GetY();       
1952         Float_t yN=0.;
1953         if (seg->GetLayer()==5) yN = yP - (7.6/1.9);
1954         else yN = yP + (7.6/1.9);
1955         seg->GetPadCxz(yP, yN, xt, zt); 
1956         
1957         if ( (TMath::Abs(xt)<hwSSD) && (TMath::Abs(zt)<hlSSD) ) {
1958           Double_t loc[3]={xt,0.,zt},trk[3]={0.,0.,0.};
1959           mT2L->MasterToLocal(loc,trk);
1960           lp[0]=trk[1];
1961           lp[1]=trk[2];        
1962           lp[4]=neg[j].GetQ(); //Q
1963           for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1964           for (Int_t ilab=0;ilab<3;ilab++) milab[ilab] = neg[j].GetLabel(ilab);   
1965           CheckLabels2(milab);
1966           milab[3]=( j << 10 ) + idet; // pos|neg|det
1967           Int_t info[3] = {0,(Int_t)neg[j].GetNd(),fNlayer[fModule]};
1968           
1969           lp[2]=7.27e-05;   // 0.0085*0.0085;  //SigmaY2
1970           lp[3]=1.33;       // 1.15*1.15;  //SigmaZ2
1971           lp[5]=0.00931;
1972           if(info[1]>1) lp[2]=6.91e-05;
1973           
1974           AliITSRecPoint * cl2;
1975           if(clusters){
1976             cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);      
1977             cl2->SetChargeRatio(1.);
1978             cl2->SetType(80);     
1979           }
1980           else{
1981             cl2 = new AliITSRecPoint(milab,lp,info);
1982             cl2->SetChargeRatio(1.);
1983             cl2->SetType(80);
1984             fDetTypeRec->AddRecPoint(*cl2);
1985           }
1986           ncl++;
1987         } // cross is within the detector
1988         //
1989         //--------------
1990         
1991       } // end loop over Pside 1D clusters
1992       
1993     } // bad Nside module
1994     
1995     //---------------------------------------------------------
1996     
1997   } // use bad channels
1998     
1999   //cout<<ncl<<" clusters for this module"<<endl;
2000
2001   delete [] cnegative;
2002   delete [] cused1;
2003   delete [] negativepair;
2004   delete [] cpositive;
2005   delete [] cused2;
2006   delete [] positivepair;
2007
2008 }
2009