]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSClusterFinderV2SSD.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
New reference data from Alexander
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderV2SSD.cxx
CommitLineData
04366a57 1/**************************************************************************
b42cfa25 2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
04366a57 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 **************************************************************************/
308b5ea4 15
16/* $Id$ */
17
04366a57 18////////////////////////////////////////////////////////////////////////////
19// Implementation of the ITS clusterer V2 class //
20// //
21// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
308b5ea4 22// Revised: Enrico Fragiacomo, enrico.fragiacomo@ts.infn.it //
04366a57 23// //
24///////////////////////////////////////////////////////////////////////////
25
308b5ea4 26#include <Riostream.h>
27
04366a57 28
29#include "AliITSClusterFinderV2SSD.h"
00a7cc50 30#include "AliITSRecPoint.h"
1f3e997f 31#include "AliITSgeomTGeo.h"
7d62fb64 32#include "AliITSDetTypeRec.h"
04366a57 33#include "AliRawReader.h"
34#include "AliITSRawStreamSSD.h"
04366a57 35#include <TClonesArray.h>
04366a57 36#include "AliITSdigitSSD.h"
a86176e3 37#include "AliITSReconstructor.h"
3a4139a2 38#include "AliITSCalibrationSSD.h"
04366a57 39
308b5ea4 40Short_t *AliITSClusterFinderV2SSD::fgPairs = 0x0;
41Int_t AliITSClusterFinderV2SSD::fgPairsSize = 0;
c157c94e 42
04366a57 43ClassImp(AliITSClusterFinderV2SSD)
44
45
e56160b8 46AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinderV2(dettyp),
1f3e997f 47fLastSSD1(AliITSgeomTGeo::GetModuleIndex(6,1,1)-1),
e56160b8 48fYpitchSSD(0.0095),
49fHwSSD(3.65),
50fHlSSD(2.00),
51fTanP(0.0275),
52fTanN(0.0075){
04366a57 53
54 //Default constructor
55
04366a57 56}
57
308b5ea4 58//______________________________________________________________________
59AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(const AliITSClusterFinderV2SSD &cf) : AliITSClusterFinderV2(cf), fLastSSD1(cf.fLastSSD1),
60fYpitchSSD(cf.fYpitchSSD),
61fHwSSD(cf.fHwSSD),
62fHlSSD(cf.fHlSSD),
63fTanP(cf.fTanP),
64fTanN(cf.fTanN)
65{
66 // Copy constructor
67
68}
69
70//______________________________________________________________________
71AliITSClusterFinderV2SSD& AliITSClusterFinderV2SSD::operator=(const AliITSClusterFinderV2SSD& cf ){
72 // Assignment operator
73
74 this->~AliITSClusterFinderV2SSD();
75 new(this) AliITSClusterFinderV2SSD(cf);
76 return *this;
77}
78
04366a57 79
80void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){
81
82 //Find clusters V2
83 SetModule(mod);
84 FindClustersSSD(fDigits);
85
86}
87
88void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) {
89 //------------------------------------------------------------
90 // Actual SSD cluster finder
91 //------------------------------------------------------------
a86176e3 92
93 static AliITSRecoParam *repa = NULL;
94 if(!repa){
95 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
96 if(!repa){
97 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParamDefault();
98 AliWarning("Using default AliITSRecoParam class");
99 }
100 }
101
3a4139a2 102 AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
103 Float_t gain=0;
104
04366a57 105 Int_t smaxall=alldigits->GetEntriesFast();
106 if (smaxall==0) return;
f7b30404 107 // TObjArray *digits = new TObjArray;
108 TObjArray digits;
04366a57 109 for (Int_t i=0;i<smaxall; i++){
110 AliITSdigitSSD *d=(AliITSdigitSSD*)alldigits->UncheckedAt(i);
3a4139a2 111
112 if(d->IsSideP()) gain = cal->GetGainP(d->GetStripNumber());
113 else gain = cal->GetGainN(d->GetStripNumber());
114
115 Float_t q=gain*d->GetSignal(); // calibration brings mip peaks around 120 (in ADC units)
116 q=cal->ADCToKeV(q); // converts the charge in KeV from ADC units
117 //Float_t q=d->GetSignal()/4.29;// temp. fix (for PID purposed - normalis. to be checked)
6490e1c5 118 d->SetSignal(Int_t(q));
3a4139a2 119
04366a57 120 if (d->GetSignal()<3) continue;
f7b30404 121 digits.AddLast(d);
04366a57 122 }
f7b30404 123 Int_t smax = digits.GetEntriesFast();
04366a57 124 if (smax==0) return;
125
126 const Int_t kMax=1000;
127 Int_t np=0, nn=0;
128 Ali1Dcluster pos[kMax], neg[kMax];
129 Float_t y=0., q=0., qmax=0.;
130 Int_t lab[4]={-2,-2,-2,-2};
131
f7b30404 132 AliITSdigitSSD *d=(AliITSdigitSSD*)digits.UncheckedAt(0);
04366a57 133 q += d->GetSignal();
134 y += d->GetCoord2()*d->GetSignal();
135 qmax=d->GetSignal();
136 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
137 Int_t curr=d->GetCoord2();
138 Int_t flag=d->GetCoord1();
139 Int_t *n=&nn;
140 Ali1Dcluster *c=neg;
141 Int_t nd=1;
142 Int_t milab[10];
143 for (Int_t ilab=0;ilab<10;ilab++){
144 milab[ilab]=-2;
145 }
146 milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2);
147
148 for (Int_t s=1; s<smax; s++) {
f7b30404 149 d=(AliITSdigitSSD*)digits.UncheckedAt(s);
04366a57 150 Int_t strip=d->GetCoord2();
151 if ((strip-curr) > 1 || flag!=d->GetCoord1()) {
152 c[*n].SetY(y/q);
153 c[*n].SetQ(q);
154 c[*n].SetNd(nd);
155 CheckLabels2(milab);
156 c[*n].SetLabels(milab);
a86176e3 157
158 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
159
160 //Split suspiciously big cluster
161 if (nd>4&&nd<25) {
162 c[*n].SetY(y/q-0.25*nd);
163 c[*n].SetQ(0.5*q);
164 (*n)++;
165 if (*n==kMax) {
166 Error("FindClustersSSD","Too many 1D clusters !");
167 return;
168 }
169 c[*n].SetY(y/q+0.25*nd);
170 c[*n].SetQ(0.5*q);
171 c[*n].SetNd(nd);
172 c[*n].SetLabels(milab);
173 }
174
175 } // unfolding is on
176
04366a57 177 (*n)++;
178 if (*n==kMax) {
179 Error("FindClustersSSD","Too many 1D clusters !");
180 return;
181 }
182 y=q=qmax=0.;
183 nd=0;
184 lab[0]=lab[1]=lab[2]=-2;
185 //
186 for (Int_t ilab=0;ilab<10;ilab++){
187 milab[ilab]=-2;
188 }
189 //
190 if (flag!=d->GetCoord1()) { n=&np; c=pos; }
191 }
192 flag=d->GetCoord1();
193 q += d->GetSignal();
194 y += d->GetCoord2()*d->GetSignal();
195 nd++;
196 if (d->GetSignal()>qmax) {
197 qmax=d->GetSignal();
198 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
199 }
200 for (Int_t ilab=0;ilab<10;ilab++) {
201 if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab)));
202 }
203 curr=strip;
204 }
205 c[*n].SetY(y/q);
206 c[*n].SetQ(q);
207 c[*n].SetNd(nd);
208 c[*n].SetLabels(lab);
a86176e3 209
210 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
211
212 //Split suspiciously big cluster
213 if (nd>4 && nd<25) {
214 c[*n].SetY(y/q-0.25*nd);
215 c[*n].SetQ(0.5*q);
216 (*n)++;
217 if (*n==kMax) {
04366a57 218 Error("FindClustersSSD","Too many 1D clusters !");
219 return;
a86176e3 220 }
221 c[*n].SetY(y/q+0.25*nd);
222 c[*n].SetQ(0.5*q);
223 c[*n].SetNd(nd);
224 c[*n].SetLabels(lab);
225 }
226 } // unfolding is on
227
04366a57 228 (*n)++;
229 if (*n==kMax) {
230 Error("FindClustersSSD","Too many 1D clusters !");
231 return;
232 }
233
234 FindClustersSSD(neg, nn, pos, np);
235}
236
237
238void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){
239
240 //------------------------------------------------------------
241 // This function creates ITS clusters from raw data
242 //------------------------------------------------------------
243 rawReader->Reset();
b42cfa25 244 /*
245 const UInt_t *evid; evid = rawReader->GetEventId();
246 cout<<"Event="<<evid[0]<<endl;
247 */
04366a57 248 AliITSRawStreamSSD inputSSD(rawReader);
8297cd8f 249 // rawReader->SelectEquipment(-1,0,15);
04366a57 250 FindClustersSSD(&inputSSD,clusters);
251
252}
253
3a4139a2 254void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStreamSSD* input,
04366a57 255 TClonesArray** clusters)
256{
257 //------------------------------------------------------------
258 // Actual SSD cluster finder for raw data
259 //------------------------------------------------------------
a86176e3 260
261 static AliITSRecoParam *repa = NULL;
262 if(!repa){
263 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
264 if(!repa){
265 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParamDefault();
266 AliWarning("Using default AliITSRecoParam class");
267 }
268 }
269
04366a57 270 Int_t nClustersSSD = 0;
271 const Int_t kMax = 1000;
272 Ali1Dcluster clusters1D[2][kMax];
273 Int_t nClusters[2] = {0, 0};
274 Int_t lab[3]={-2,-2,-2};
275 Float_t q = 0.;
276 Float_t y = 0.;
277 Int_t nDigits = 0;
3a4139a2 278 Float_t gain=0;
308b5ea4 279 Float_t noise=0.;
5ffe058f 280 // Float_t pedestal=0.;
308b5ea4 281 Float_t oldnoise=0.;
b4704be3 282 AliITSCalibrationSSD* cal=NULL;
308b5ea4 283
308b5ea4 284 Int_t matrix[12][1536];
d5ad0697 285 Int_t iddl=-1;
286 Int_t iad=-1;
287 Int_t oddl = -1;
288 Int_t oad = -1;
289 Int_t oadc = -1;
290 Int_t ostrip = -1;
308b5ea4 291 Int_t osignal = 65535;
308b5ea4 292 Int_t n=0;
293 Bool_t next=0;
04366a57 294
295 // read raw data input stream
296 while (kTRUE) {
308b5ea4 297
298 // reset signal matrix
299 for(Int_t i=0; i<12; i++) { for(Int_t j=0; j<1536; j++) { matrix[i][j] = 65535;} }
b42cfa25 300
308b5ea4 301 if(osignal!=65535) {
302 n++;
303 matrix[oadc][ostrip] = osignal; // recover data from previous occurence of input->Next()
304 }
b42cfa25 305
308b5ea4 306 // buffer data for ddl=iddl and ad=iad
307 while(kTRUE) {
bc4dd89a 308
308b5ea4 309 next = input->Next();
d5ad0697 310 if((!next)&&(input->flag)) continue;
311 Int_t ddl=input->GetDDL();
312 Int_t ad=input->GetAD();
313 Int_t adc = input->GetADC(); adc = (adc<6)? adc : adc - 2;
314 Int_t strip = input->GetStrip();
308b5ea4 315 if(input->GetSideFlag()) strip=1535-strip;
d5ad0697 316 Int_t signal = input->GetSignal();
308b5ea4 317 //cout<<ddl<<" "<<ad<<" "<<adc<<" "<<strip<<" "<<signal<<endl;
b42cfa25 318
308b5ea4 319 if((ddl==iddl)&&(ad==iad)) {n++; matrix[adc][strip] = signal;}
320 else {oddl=iddl; oad=iad; oadc = adc; ostrip = strip; osignal=signal; iddl=ddl; iad=ad; break;}
321
b42cfa25 322 if(!next) {oddl=iddl; oad=iad; oadc = adc; ostrip = strip; osignal=signal; iddl=ddl; iad=ad; break;}
323 //break;
308b5ea4 324 }
b42cfa25 325
326 // No SSD data
d5ad0697 327 if(!next && oddl<0) break;
b42cfa25 328
308b5ea4 329 if(n==0) continue; // first occurence
330 n=0; osignal=0;
b42cfa25 331
308b5ea4 332 // fill 1Dclusters
333 for(Int_t iadc=0; iadc<12; iadc++) { // loop over ADC index for ddl=oddl and ad=oad
b42cfa25 334
308b5ea4 335 Int_t iimod = (oad - 1) * 12 + iadc;
d5ad0697 336 Int_t iModule = AliITSRawStreamSSD::GetModuleNumber(oddl,iimod);
308b5ea4 337 if(iModule==-1) continue;
b42cfa25 338 // cout<<"ddl="<<oddl<<" ad"<<oad<<" module="<<iModule<<endl;
308b5ea4 339 cal = (AliITSCalibrationSSD*)GetResp(iModule);
b42cfa25 340
308b5ea4 341 Bool_t first = 0;
b42cfa25 342
bc4dd89a 343 /*
308b5ea4 344 for(Int_t istrip=0; istrip<768; istrip++) { // P-side
d5ad0697 345 Int_t signal = matrix[iadc][istrip];
b42cfa25 346 pedestal = cal->GetPedestalP(istrip);
e7f0e76b 347 matrix[iadc][istrip]=signal-(Int_t)pedestal;
bc4dd89a 348 }
349 */
350
351 /*
b42cfa25 352 Float_t cmode=0;
353 for(Int_t l=0; l<6; l++) {
354 cmode=0;
355 for(Int_t n=20; n<108; n++) cmode+=matrix[iadc][l*128+n];
356 cmode/=88.;
e7f0e76b 357 for(Int_t n=0; n<128; n++) matrix[iadc][l*128+n]-=(Int_t)cmode;
b42cfa25 358
359 }
bc4dd89a 360 */
361
b42cfa25 362 for(Int_t istrip=0; istrip<768; istrip++) { // P-side
363
364 Int_t signal = TMath::Abs(matrix[iadc][istrip]);
365
308b5ea4 366 oldnoise = noise;
b42cfa25 367 noise = cal->GetNoiseP(istrip); if(noise<1.) signal = 65535;
368 if(signal<5*noise) signal = 65535; // in case ZS was not done in hw do it now
369 if( (signal<30.) || (istrip<10) || (istrip>758) ) signal=65535;
308b5ea4 370
371 if (signal!=65535) {
372 gain = cal->GetGainP(istrip);
373 signal = (Int_t) ( signal * gain ); // signal is corrected for gain
308b5ea4 374 signal = (Int_t) cal->ADCToKeV( signal ); // signal is converted in KeV
b42cfa25 375
308b5ea4 376 q += signal; // add digit to current cluster
377 y += istrip * signal;
378 nDigits++;
379 first=1;
380 }
b42cfa25 381
308b5ea4 382 else if(first) {
b42cfa25 383
308b5ea4 384 if ( ((nDigits==1)&&(q>5*oldnoise)) || (nDigits>1) ) {
b42cfa25 385
308b5ea4 386 Ali1Dcluster& cluster = clusters1D[0][nClusters[0]++];
387 cluster.SetY(y/q);
388 cluster.SetQ(q);
389 cluster.SetNd(nDigits);
390 cluster.SetLabels(lab);
391
a86176e3 392 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
393
394 //Split suspiciously big cluster
395 if (nDigits > 4&&nDigits < 25) {
396 cluster.SetY(y/q - 0.25*nDigits);
397 cluster.SetQ(0.5*q);
398 if (nClusters[0] == kMax) {
399 Error("FindClustersSSD", "Too many 1D clusters !");
400 return;
401 }
402 Ali1Dcluster& cluster2 = clusters1D[0][nClusters[0]++];
403 cluster2.SetY(y/q + 0.25*nDigits);
404 cluster2.SetQ(0.5*q);
405 cluster2.SetNd(nDigits);
406 cluster2.SetLabels(lab);
308b5ea4 407 }
a86176e3 408 } // unfolding is on
308b5ea4 409 }
410
411 y = q = 0.;
412 nDigits = 0;
413 first=0;
04366a57 414 }
a86176e3 415
308b5ea4 416 } // loop over strip on P-side
a86176e3 417
308b5ea4 418 // if last strip does have signal
419 if(first) {
420
a86176e3 421 if ( ((nDigits==1)&&(q>5*oldnoise)) || (nDigits>1) ) {
422
423 Ali1Dcluster& cluster = clusters1D[0][nClusters[0]++];
424 cluster.SetY(y/q);
425 cluster.SetQ(q);
426 cluster.SetNd(nDigits);
427 cluster.SetLabels(lab);
428
429 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
308b5ea4 430
431 //Split suspiciously big cluster
432 if (nDigits > 4&&nDigits < 25) {
433 cluster.SetY(y/q - 0.25*nDigits);
434 cluster.SetQ(0.5*q);
435 if (nClusters[0] == kMax) {
436 Error("FindClustersSSD", "Too many 1D clusters !");
437 return;
438 }
439 Ali1Dcluster& cluster2 = clusters1D[0][nClusters[0]++];
440 cluster2.SetY(y/q + 0.25*nDigits);
441 cluster2.SetQ(0.5*q);
442 cluster2.SetNd(nDigits);
443 cluster2.SetLabels(lab);
444 }
a86176e3 445 } // unfolding is on
446
447 }
448 y = q = 0.;
449 nDigits = 0;
450 first=0;
04366a57 451 }
308b5ea4 452
bc4dd89a 453 /*
b42cfa25 454 for(Int_t istrip=768; istrip<1536; istrip++) { // P-side
455 Int_t signal = matrix[iadc][istrip];
456 pedestal = cal->GetPedestalN(1535-istrip);
e7f0e76b 457 matrix[iadc][istrip]=signal-(Int_t)pedestal;
b42cfa25 458 }
bc4dd89a 459 */
b42cfa25 460
bc4dd89a 461 /*
b42cfa25 462 for(Int_t l=6; l<12; l++) {
463 Float_t cmode=0;
464 for(Int_t n=20; n<108; n++) cmode+=matrix[iadc][l*128+n];
465 cmode/=88.;
e7f0e76b 466 for(Int_t n=0; n<128; n++) matrix[iadc][l*128+n]-=(Int_t)cmode;
b42cfa25 467 }
bc4dd89a 468 */
b42cfa25 469
308b5ea4 470 oldnoise = 0.;
471 noise = 0.;
472 for(Int_t istrip=768; istrip<1536; istrip++) { // N-side
473
b42cfa25 474 Int_t signal = TMath::Abs(matrix[iadc][istrip]);
475 //cout<<"####"<<" "<<oddl<<" "<<oad<<" "<<iadc<<" "<<istrip<<" "<<signal<<endl;
476
d5ad0697 477 Int_t strip = 1535-istrip;
b42cfa25 478
308b5ea4 479 oldnoise = noise;
b42cfa25 480 noise = cal->GetNoiseN(strip); if(noise<1.) signal=65535;
481 if(signal<5*noise) signal = 65535; // in case ZS was not done in hw do it now
482 if( (signal<30.) || (istrip<778) || (istrip>1526) ) signal=65535;
308b5ea4 483
484 if (signal!=65535) {
b42cfa25 485 // cout<<"ddl="<<oddl<<" ad"<<oad<<" module="<<iModule<<" strip= "<<istrip<<
486 // " sig="<<signal<<" "<<cal->GetPedestalN(strip)<<endl;
308b5ea4 487 gain = cal->GetGainN(strip);
488 signal = (Int_t) ( signal * gain); // signal is corrected for gain
489 signal = (Int_t) cal->ADCToKeV( signal ); // signal is converted in KeV
490
491 // add digit to current cluster
492 q += signal;
493 y += strip * signal;
494 nDigits++;
495 first=1;
496 }
497
498 else if(first) {
a86176e3 499
308b5ea4 500 if ( ((nDigits==1)&&(q>5*oldnoise)) || (nDigits>1) ) {
501
502 Ali1Dcluster& cluster = clusters1D[1][nClusters[1]++];
503 cluster.SetY(y/q);
504 cluster.SetQ(q);
505 cluster.SetNd(nDigits);
506 cluster.SetLabels(lab);
507
a86176e3 508 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
509
510 //Split suspiciously big cluster
511 if (nDigits > 4&&nDigits < 25) {
512 cluster.SetY(y/q - 0.25*nDigits);
513 cluster.SetQ(0.5*q);
514 if (nClusters[1] == kMax) {
515 Error("FindClustersSSD", "Too many 1D clusters !");
516 return;
517 }
518 Ali1Dcluster& cluster2 = clusters1D[1][nClusters[1]++];
519 cluster2.SetY(y/q + 0.25*nDigits);
520 cluster2.SetQ(0.5*q);
521 cluster2.SetNd(nDigits);
522 cluster2.SetLabels(lab);
523 }
524 } // unfolding is on
525 }
526
308b5ea4 527 y = q = 0.;
528 nDigits = 0;
529 first=0;
530 }
531
532 } // loop over strips on N-side
04366a57 533
308b5ea4 534 if(first) {
535
a86176e3 536 if ( ((nDigits==1)&&(q>5*oldnoise)) || (nDigits>1) ) {
537
538 Ali1Dcluster& cluster = clusters1D[1][nClusters[1]++];
539 cluster.SetY(y/q);
540 cluster.SetQ(q);
541 cluster.SetNd(nDigits);
542 cluster.SetLabels(lab);
543
544 if(repa->GetUseUnfoldingInClusterFinderSSD()==kTRUE) {
308b5ea4 545
546 //Split suspiciously big cluster
547 if (nDigits > 4&&nDigits < 25) {
548 cluster.SetY(y/q - 0.25*nDigits);
549 cluster.SetQ(0.5*q);
550 if (nClusters[1] == kMax) {
551 Error("FindClustersSSD", "Too many 1D clusters !");
552 return;
553 }
554 Ali1Dcluster& cluster2 = clusters1D[1][nClusters[1]++];
555 cluster2.SetY(y/q + 0.25*nDigits);
556 cluster2.SetQ(0.5*q);
557 cluster2.SetNd(nDigits);
558 cluster2.SetLabels(lab);
559 }
a86176e3 560 } // unfolding is on
561 }
562
563 y = q = 0.;
564 nDigits = 0;
565 first=0;
308b5ea4 566 }
567
568 // create recpoints
569 if((nClusters[0])&&(nClusters[1])) {
a86176e3 570
308b5ea4 571 //cout<<"creating recpoint for module="<<iModule<<" "<<nClusters[0]<<" "<<nClusters[1]<<endl;
00a7cc50 572 clusters[iModule] = new TClonesArray("AliITSRecPoint");
04366a57 573 fModule = iModule;
308b5ea4 574 // fModule = 500;
04366a57 575 FindClustersSSD(&clusters1D[0][0], nClusters[0],
576 &clusters1D[1][0], nClusters[1], clusters[iModule]);
577 Int_t nClusters = clusters[iModule]->GetEntriesFast();
578 nClustersSSD += nClusters;
579 }
580
04366a57 581 nClusters[0] = nClusters[1] = 0;
582 y = q = 0.;
583 nDigits = 0;
3a4139a2 584
308b5ea4 585 } // loop over adc
04366a57 586
d5ad0697 587 if(!next) break;
04366a57 588 }
308b5ea4 589
04366a57 590 Info("FindClustersSSD", "found clusters in ITS SSD: %d", nClustersSSD);
591}
592
593void AliITSClusterFinderV2SSD::
594FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
595 Ali1Dcluster* pos, Int_t np,
596 TClonesArray *clusters) {
597 //------------------------------------------------------------
598 // Actual SSD cluster finder
599 //------------------------------------------------------------
b4704be3 600
308b5ea4 601 // Float_t xyz[3];
602
b4704be3 603 const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule);
604
04366a57 605 TClonesArray &cl=*clusters;
606 //
607 Float_t tanp=fTanP, tann=fTanN;
608 if (fModule>fLastSSD1) {tann=fTanP; tanp=fTanN;}
609 Int_t idet=fNdet[fModule];
610 Int_t ncl=0;
611 //
612 Int_t negativepair[30000];
613 Int_t cnegative[3000];
614 Int_t cused1[3000];
615 Int_t positivepair[30000];
616 Int_t cpositive[3000];
617 Int_t cused2[3000];
618 for (Int_t i=0;i<3000;i++) {cnegative[i]=0; cused1[i]=0;}
619 for (Int_t i=0;i<3000;i++) {cpositive[i]=0; cused2[i]=0;}
78538bfe 620 for (Int_t i=0;i<30000;i++) {negativepair[i]=0; positivepair[i]=0;}
c157c94e 621
308b5ea4 622 if ((np*nn) > fgPairsSize) {
623 if (fgPairs) delete [] fgPairs;
624 fgPairsSize = 4*np*nn;
625 fgPairs = new Short_t[fgPairsSize];
c157c94e 626 }
308b5ea4 627 memset(fgPairs,0,sizeof(Short_t)*np*nn);
628
04366a57 629 //
630 // find available pairs
631 //
632 for (Int_t i=0; i<np; i++) {
633 Float_t yp=pos[i].GetY()*fYpitchSSD;
634 if (pos[i].GetQ()<3) continue;
635 for (Int_t j=0; j<nn; j++) {
636 if (neg[j].GetQ()<3) continue;
637 Float_t yn=neg[j].GetY()*fYpitchSSD;
638 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
639 Float_t yt=yn + tann*zt;
640 zt-=fHlSSD; yt-=fHwSSD;
641 if (TMath::Abs(yt)<fHwSSD+0.01)
642 if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
643 negativepair[i*10+cnegative[i]] =j; //index
644 positivepair[j*10+cpositive[j]] =i;
645 cnegative[i]++; //counters
646 cpositive[j]++;
308b5ea4 647 fgPairs[i*nn+j]=100;
04366a57 648 }
649 }
650 }
308b5ea4 651
652 //
653 // try to recover points out of but close to the module boundaries
04366a57 654 //
655 for (Int_t i=0; i<np; i++) {
656 Float_t yp=pos[i].GetY()*fYpitchSSD;
657 if (pos[i].GetQ()<3) continue;
658 for (Int_t j=0; j<nn; j++) {
659 if (neg[j].GetQ()<3) continue;
308b5ea4 660 // if both 1Dclusters have an other cross continue
04366a57 661 if (cpositive[j]&&cnegative[i]) continue;
662 Float_t yn=neg[j].GetY()*fYpitchSSD;
663 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
664 Float_t yt=yn + tann*zt;
665 zt-=fHlSSD; yt-=fHwSSD;
666 if (TMath::Abs(yt)<fHwSSD+0.1)
667 if (TMath::Abs(zt)<fHlSSD+0.15) {
308b5ea4 668 // tag 1Dcluster (eventually will produce low quality recpoint)
04366a57 669 if (cnegative[i]==0) pos[i].SetNd(100); // not available pair
670 if (cpositive[j]==0) neg[j].SetNd(100); // not available pair
671 negativepair[i*10+cnegative[i]] =j; //index
672 positivepair[j*10+cpositive[j]] =i;
673 cnegative[i]++; //counters
674 cpositive[j]++;
308b5ea4 675 fgPairs[i*nn+j]=100;
04366a57 676 }
677 }
678 }
308b5ea4 679
04366a57 680 //
681 Float_t lp[5];
682 Int_t milab[10];
683 Double_t ratio;
684
b42cfa25 685
a86176e3 686 static AliITSRecoParam *repa = NULL;
687 if(!repa){
688 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
689 if(!repa){
690 repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParamDefault();
691 AliWarning("Using default AliITSRecoParam class");
692 }
693 }
b42cfa25 694
a86176e3 695 if(repa->GetUseChargeMatchingInClusterFinderSSD()==kTRUE) {
308b5ea4 696
697
308b5ea4 698
308b5ea4 699
a86176e3 700 //
701 // sign gold tracks
702 //
703 for (Int_t ip=0;ip<np;ip++){
704 Float_t ybest=1000,zbest=1000,qbest=0;
705 //
706 // select gold clusters
707 if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){
708 Float_t yp=pos[ip].GetY()*fYpitchSSD;
709 Int_t j = negativepair[10*ip];
710 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
711 //
712 Float_t yn=neg[j].GetY()*fYpitchSSD;
713 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
714 Float_t yt=yn + tann*zt;
715 zt-=fHlSSD; yt-=fHwSSD;
716 ybest=yt; zbest=zt;
717 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
718
719 //cout<<yt<<" "<<zt<<" "<<qbest<<endl;
720
721 {
722 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
723 mT2L->MasterToLocal(loc,trk);
724 lp[0]=trk[1];
725 lp[1]=trk[2];
00a7cc50 726 }
a86176e3 727 lp[2]=0.0025*0.0025; //SigmaY2
728 lp[3]=0.110*0.110; //SigmaZ2
00a7cc50 729
a86176e3 730 lp[4]=qbest; //Q
731 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
732 for (Int_t ilab=0;ilab<3;ilab++){
733 milab[ilab] = pos[ip].GetLabel(ilab);
734 milab[ilab+3] = neg[j].GetLabel(ilab);
00a7cc50 735 }
04366a57 736 //
a86176e3 737 CheckLabels2(milab);
738 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
739 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
740 AliITSRecPoint * cl2;
741
742 if(clusters){ // Note clusters != 0 when method is called for rawdata
04366a57 743
a86176e3 744
745 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
746
747 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
748
749 cl2->SetChargeRatio(ratio);
750 cl2->SetType(1);
751 fgPairs[ip*nn+j]=1;
752 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
753 cl2->SetType(2);
754 fgPairs[ip*nn+j]=2;
04366a57 755 }
a86176e3 756 cused1[ip]++;
757 cused2[j]++;
758
759 }
760 else{ // Note clusters == 0 when method is called for digits
761
762 cl2 = new AliITSRecPoint(milab,lp,info);
763
764 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
765
766 cl2->SetChargeRatio(ratio);
767 cl2->SetType(1);
768 fgPairs[ip*nn+j]=1;
769 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
770 cl2->SetType(2);
771 fgPairs[ip*nn+j]=2;
00a7cc50 772 }
a86176e3 773 cused1[ip]++;
774 cused2[j]++;
775 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" gold"<<endl;
776 fDetTypeRec->AddRecPoint(*cl2);
777 }
778 ncl++;
779 }
780 }
781
782 for (Int_t ip=0;ip<np;ip++){
783 Float_t ybest=1000,zbest=1000,qbest=0;
784 //
785 //
786 // select "silber" cluster
787 if ( cnegative[ip]==1 && cpositive[negativepair[10*ip]]==2){
788 Int_t in = negativepair[10*ip];
789 Int_t ip2 = positivepair[10*in];
790 if (ip2==ip) ip2 = positivepair[10*in+1];
791 Float_t pcharge = pos[ip].GetQ()+pos[ip2].GetQ();
792 if (TMath::Abs(pcharge-neg[in].GetQ())<10){
793 //
794 // add first pair
795 if (fgPairs[ip*nn+in]==100){ //
796 Float_t yp=pos[ip].GetY()*fYpitchSSD;
797 Float_t yn=neg[in].GetY()*fYpitchSSD;
798 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
799 Float_t yt=yn + tann*zt;
800 zt-=fHlSSD; yt-=fHwSSD;
801 ybest =yt; zbest=zt;
802 qbest =pos[ip].GetQ();
803 {
804 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
805 mT2L->MasterToLocal(loc,trk);
806 lp[0]=trk[1];
807 lp[1]=trk[2];
00a7cc50 808 }
a86176e3 809 lp[2]=0.0025*0.0025; //SigmaY2
810 lp[3]=0.110*0.110; //SigmaZ2
811
812 lp[4]=qbest; //Q
813 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
814 for (Int_t ilab=0;ilab<3;ilab++){
815 milab[ilab] = pos[ip].GetLabel(ilab);
816 milab[ilab+3] = neg[in].GetLabel(ilab);
817 }
818 //
819 CheckLabels2(milab);
820 ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ());
821 milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det
822 Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
823
824 AliITSRecPoint * cl2;
825 if(clusters){
826
827 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
828
829 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
830
831 cl2->SetChargeRatio(ratio);
832 cl2->SetType(5);
833 fgPairs[ip*nn+in] = 5;
834 if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
835 cl2->SetType(6);
836 fgPairs[ip*nn+in] = 6;
837 }
838 }
839 else{
840 cl2 = new AliITSRecPoint(milab,lp,info);
841 cl2->SetChargeRatio(ratio);
842 cl2->SetType(5);
843 fgPairs[ip*nn+in] = 5;
844 if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
845 cl2->SetType(6);
846 fgPairs[ip*nn+in] = 6;
847 }
848 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" silver1"<<endl;
849
850 fDetTypeRec->AddRecPoint(*cl2);
851 }
852 ncl++;
04366a57 853 }
04366a57 854
04366a57 855 //
a86176e3 856 // add second pair
00a7cc50 857
a86176e3 858 // if (!(cused1[ip2] || cused2[in])){ //
859 if (fgPairs[ip2*nn+in]==100){
860 Float_t yp=pos[ip2].GetY()*fYpitchSSD;
861 Float_t yn=neg[in].GetY()*fYpitchSSD;
862 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
863 Float_t yt=yn + tann*zt;
864 zt-=fHlSSD; yt-=fHwSSD;
865 ybest =yt; zbest=zt;
866 qbest =pos[ip2].GetQ();
867 {
868 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
869 mT2L->MasterToLocal(loc,trk);
870 lp[0]=trk[1];
871 lp[1]=trk[2];
00a7cc50 872 }
a86176e3 873 lp[2]=0.0025*0.0025; //SigmaY2
874 lp[3]=0.110*0.110; //SigmaZ2
875
876 lp[4]=qbest; //Q
877 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
878 for (Int_t ilab=0;ilab<3;ilab++){
879 milab[ilab] = pos[ip2].GetLabel(ilab);
880 milab[ilab+3] = neg[in].GetLabel(ilab);
00a7cc50 881 }
a86176e3 882 //
883 CheckLabels2(milab);
884 ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ());
885 milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det
886 Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
308b5ea4 887
a86176e3 888 AliITSRecPoint * cl2;
889 if(clusters){
890 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
891
892 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
893
894 cl2->SetChargeRatio(ratio);
895 cl2->SetType(5);
896 fgPairs[ip2*nn+in] =5;
897 if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
898 cl2->SetType(6);
899 fgPairs[ip2*nn+in] =6;
900 }
901 }
902 else{
903 cl2 = new AliITSRecPoint(milab,lp,info);
904 cl2->SetChargeRatio(ratio);
905 cl2->SetType(5);
906 fgPairs[ip2*nn+in] =5;
907 if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
908 cl2->SetType(6);
909 fgPairs[ip2*nn+in] =6;
910 }
911
912 // cout<<"AliITSClusterFinderV2SSD "<<fModule<<" silver2"<<endl;
913 fDetTypeRec->AddRecPoint(*cl2);
914 }
915 ncl++;
916 }
917 cused1[ip]++;
918 cused1[ip2]++;
919 cused2[in]++;
920 }
921 }
922 }
923
924 } // use charge matching
925
926
04366a57 927 //
928 for (Int_t jn=0;jn<nn;jn++){
929 if (cused2[jn]) continue;
930 Float_t ybest=1000,zbest=1000,qbest=0;
931 // select "silber" cluster
932 if ( cpositive[jn]==1 && cnegative[positivepair[10*jn]]==2){
933 Int_t ip = positivepair[10*jn];
934 Int_t jn2 = negativepair[10*ip];
935 if (jn2==jn) jn2 = negativepair[10*ip+1];
936 Float_t pcharge = neg[jn].GetQ()+neg[jn2].GetQ();
937 //
938 if (TMath::Abs(pcharge-pos[ip].GetQ())<10){
939 //
940 // add first pair
941 // if (!(cused1[ip]||cused2[jn])){
308b5ea4 942 if (fgPairs[ip*nn+jn]==100){
04366a57 943 Float_t yn=neg[jn].GetY()*fYpitchSSD;
944 Float_t yp=pos[ip].GetY()*fYpitchSSD;
945 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
946 Float_t yt=yn + tann*zt;
947 zt-=fHlSSD; yt-=fHwSSD;
948 ybest =yt; zbest=zt;
949 qbest =neg[jn].GetQ();
b4704be3 950 {
951 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
952 mT2L->MasterToLocal(loc,trk);
953 lp[0]=trk[1];
954 lp[1]=trk[2];
955 }
04366a57 956 lp[2]=0.0025*0.0025; //SigmaY2
957 lp[3]=0.110*0.110; //SigmaZ2
958
959 lp[4]=qbest; //Q
960 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
961 for (Int_t ilab=0;ilab<3;ilab++){
962 milab[ilab] = pos[ip].GetLabel(ilab);
963 milab[ilab+3] = neg[jn].GetLabel(ilab);
964 }
965 //
966 CheckLabels2(milab);
967 ratio = (pos[ip].GetQ()-neg[jn].GetQ())/(pos[ip].GetQ()+neg[jn].GetQ());
968 milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det
969 Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]};
970
00a7cc50 971 AliITSRecPoint * cl2;
972 if(clusters){
75fb37cc 973 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 974
975 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
976
00a7cc50 977 cl2->SetChargeRatio(ratio);
978 cl2->SetType(7);
308b5ea4 979 fgPairs[ip*nn+jn] =7;
00a7cc50 980 if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
981 cl2->SetType(8);
308b5ea4 982 fgPairs[ip*nn+jn]=8;
00a7cc50 983 }
984
985 }
04366a57 986 else{
75fb37cc 987 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 988 cl2->SetChargeRatio(ratio);
989 cl2->SetType(7);
308b5ea4 990 fgPairs[ip*nn+jn] =7;
00a7cc50 991 if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
992 cl2->SetType(8);
308b5ea4 993 fgPairs[ip*nn+jn]=8;
00a7cc50 994 }
308b5ea4 995 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" silverN1"<<endl;
00a7cc50 996
997 fDetTypeRec->AddRecPoint(*cl2);
04366a57 998 }
999 ncl++;
04366a57 1000 }
1001 //
1002 // add second pair
1003 // if (!(cused1[ip]||cused2[jn2])){
308b5ea4 1004 if (fgPairs[ip*nn+jn2]==100){
04366a57 1005 Float_t yn=neg[jn2].GetY()*fYpitchSSD;
1006 Double_t yp=pos[ip].GetY()*fYpitchSSD;
1007 Double_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
1008 Double_t yt=yn + tann*zt;
1009 zt-=fHlSSD; yt-=fHwSSD;
1010 ybest =yt; zbest=zt;
1011 qbest =neg[jn2].GetQ();
b4704be3 1012 {
1013 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
1014 mT2L->MasterToLocal(loc,trk);
1015 lp[0]=trk[1];
1016 lp[1]=trk[2];
1017 }
04366a57 1018 lp[2]=0.0025*0.0025; //SigmaY2
1019 lp[3]=0.110*0.110; //SigmaZ2
1020
1021 lp[4]=qbest; //Q
1022 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1023 for (Int_t ilab=0;ilab<3;ilab++){
1024 milab[ilab] = pos[ip].GetLabel(ilab);
1025 milab[ilab+3] = neg[jn2].GetLabel(ilab);
1026 }
1027 //
1028 CheckLabels2(milab);
1029 ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ());
1030 milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det
1031 Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]};
00a7cc50 1032 AliITSRecPoint * cl2;
1033 if(clusters){
75fb37cc 1034 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1035
1036 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
1037
00a7cc50 1038 cl2->SetChargeRatio(ratio);
308b5ea4 1039 fgPairs[ip*nn+jn2]=7;
00a7cc50 1040 cl2->SetType(7);
1041 if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
1042 cl2->SetType(8);
308b5ea4 1043 fgPairs[ip*nn+jn2]=8;
00a7cc50 1044 }
1045
1046 }
04366a57 1047 else{
75fb37cc 1048 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 1049 cl2->SetChargeRatio(ratio);
308b5ea4 1050 fgPairs[ip*nn+jn2]=7;
00a7cc50 1051 cl2->SetType(7);
1052 if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
1053 cl2->SetType(8);
308b5ea4 1054 fgPairs[ip*nn+jn2]=8;
00a7cc50 1055 }
308b5ea4 1056 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" silverN2"<<endl;
1057
00a7cc50 1058 fDetTypeRec->AddRecPoint(*cl2);
04366a57 1059 }
1060
1061 ncl++;
04366a57 1062 }
1063 cused1[ip]++;
1064 cused2[jn]++;
1065 cused2[jn2]++;
1066 }
1067 }
1068 }
1069
1070 for (Int_t ip=0;ip<np;ip++){
1071 Float_t ybest=1000,zbest=1000,qbest=0;
1072 //
1073 // 2x2 clusters
1074 //
1075 if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){
1076 Float_t minchargediff =4.;
1077 Int_t j=-1;
1078 for (Int_t di=0;di<cnegative[ip];di++){
1079 Int_t jc = negativepair[ip*10+di];
1080 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1081 if (TMath::Abs(chargedif)<minchargediff){
1082 j =jc;
1083 minchargediff = TMath::Abs(chargedif);
1084 }
1085 }
1086 if (j<0) continue; // not proper cluster
308b5ea4 1087
04366a57 1088 Int_t count =0;
1089 for (Int_t di=0;di<cnegative[ip];di++){
1090 Int_t jc = negativepair[ip*10+di];
1091 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
1092 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1093 }
1094 if (count>1) continue; // more than one "proper" cluster for positive
1095 //
1096 count =0;
1097 for (Int_t dj=0;dj<cpositive[j];dj++){
1098 Int_t ic = positivepair[j*10+dj];
1099 Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
1100 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
1101 }
1102 if (count>1) continue; // more than one "proper" cluster for negative
1103
1104 Int_t jp = 0;
1105
1106 count =0;
1107 for (Int_t dj=0;dj<cnegative[jp];dj++){
1108 Int_t ic = positivepair[jp*10+dj];
1109 Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
1110 if (TMath::Abs(chargedif)<minchargediff+4.) count++;
1111 }
1112 if (count>1) continue;
308b5ea4 1113 if (fgPairs[ip*nn+j]<100) continue;
04366a57 1114 //
1115 //almost gold clusters
1116 Float_t yp=pos[ip].GetY()*fYpitchSSD;
1117 Float_t yn=neg[j].GetY()*fYpitchSSD;
1118 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
1119 Float_t yt=yn + tann*zt;
1120 zt-=fHlSSD; yt-=fHwSSD;
1121 ybest=yt; zbest=zt;
1122 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
b4704be3 1123 {
1124 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
1125 mT2L->MasterToLocal(loc,trk);
1126 lp[0]=trk[1];
1127 lp[1]=trk[2];
1128 }
04366a57 1129 lp[2]=0.0025*0.0025; //SigmaY2
1130 lp[3]=0.110*0.110; //SigmaZ2
1131 lp[4]=qbest; //Q
1132 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1133 for (Int_t ilab=0;ilab<3;ilab++){
1134 milab[ilab] = pos[ip].GetLabel(ilab);
1135 milab[ilab+3] = neg[j].GetLabel(ilab);
1136 }
1137 //
1138 CheckLabels2(milab);
1139 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
1140 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
1141 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
00a7cc50 1142 AliITSRecPoint * cl2;
1143 if(clusters){
75fb37cc 1144 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1145
1146 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
1147
00a7cc50 1148 cl2->SetChargeRatio(ratio);
1149 cl2->SetType(10);
308b5ea4 1150 fgPairs[ip*nn+j]=10;
00a7cc50 1151 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1152 cl2->SetType(11);
308b5ea4 1153 fgPairs[ip*nn+j]=11;
00a7cc50 1154 }
1155 cused1[ip]++;
1156 cused2[j]++;
04366a57 1157 }
00a7cc50 1158 else{
75fb37cc 1159 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 1160 cl2->SetChargeRatio(ratio);
1161 cl2->SetType(10);
308b5ea4 1162 fgPairs[ip*nn+j]=10;
00a7cc50 1163 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
1164 cl2->SetType(11);
308b5ea4 1165 fgPairs[ip*nn+j]=11;
00a7cc50 1166 }
1167 cused1[ip]++;
1168 cused2[j]++;
1169
308b5ea4 1170 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" 2x2"<<endl;
1171
00a7cc50 1172 fDetTypeRec->AddRecPoint(*cl2);
1173 }
04366a57 1174 ncl++;
04366a57 1175 }
1176
1177 }
1178
1179 //
1180 for (Int_t i=0; i<np; i++) {
1181 Float_t ybest=1000,zbest=1000,qbest=0;
1182 Float_t yp=pos[i].GetY()*fYpitchSSD;
1183 if (pos[i].GetQ()<3) continue;
1184 for (Int_t j=0; j<nn; j++) {
1185 // for (Int_t di = 0;di<cpositive[i];di++){
1186 // Int_t j = negativepair[10*i+di];
1187 if (neg[j].GetQ()<3) continue;
1188 if (cused2[j]||cused1[i]) continue;
308b5ea4 1189 if (fgPairs[i*nn+j]>0 &&fgPairs[i*nn+j]<100) continue;
04366a57 1190 ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ());
1191 Float_t yn=neg[j].GetY()*fYpitchSSD;
1192 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
1193 Float_t yt=yn + tann*zt;
1194 zt-=fHlSSD; yt-=fHwSSD;
1195 if (TMath::Abs(yt)<fHwSSD+0.01)
1196 if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
1197 ybest=yt; zbest=zt;
1198 qbest=0.5*(pos[i].GetQ()+neg[j].GetQ());
b4704be3 1199 {
1200 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
1201 mT2L->MasterToLocal(loc,trk);
1202 lp[0]=trk[1];
1203 lp[1]=trk[2];
1204 }
04366a57 1205 lp[2]=0.0025*0.0025; //SigmaY2
1206 lp[3]=0.110*0.110; //SigmaZ2
1207
1208 lp[4]=qbest; //Q
1209 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
1210 for (Int_t ilab=0;ilab<3;ilab++){
1211 milab[ilab] = pos[i].GetLabel(ilab);
1212 milab[ilab+3] = neg[j].GetLabel(ilab);
1213 }
1214 //
1215 CheckLabels2(milab);
1216 milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det
1217 Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
00a7cc50 1218 AliITSRecPoint * cl2;
1219 if(clusters){
75fb37cc 1220 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
308b5ea4 1221
1222 // cl2-> GetGlobalXYZ(xyz); cout<<"rec "<<xyz[0]<<" "<<xyz[1]<<" "<<xyz[2]<<endl;
1223
00a7cc50 1224 cl2->SetChargeRatio(ratio);
1225 cl2->SetType(100+cpositive[j]+cnegative[i]);
1226 }
04366a57 1227 else{
75fb37cc 1228 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 1229 cl2->SetChargeRatio(ratio);
1230 cl2->SetType(100+cpositive[j]+cnegative[i]);
308b5ea4 1231
1232 //cout<<"AliITSClusterFinderV2SSD "<<fModule<<" other"<<endl;
1233
00a7cc50 1234 fDetTypeRec->AddRecPoint(*cl2);
04366a57 1235 }
1236 ncl++;
04366a57 1237 //cl2->SetType(0);
1238 /*
308b5ea4 1239 if (fgPairs[i*nn+j]<100){
1240 printf("problem:- %d\n", fgPairs[i*nn+j]);
04366a57 1241 }
1242 if (cnegative[i]<2&&cpositive[j]<2){
308b5ea4 1243 printf("problem:- %d\n", fgPairs[i*nn+j]);
04366a57 1244 }
1245 */
1246 }
1247 }
1248 }
1249
04366a57 1250}
1251