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