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