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