]> 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
Part of Reset() moved to base class.
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderV2SSD.cxx
CommitLineData
04366a57 1/**************************************************************************
2 * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15////////////////////////////////////////////////////////////////////////////
16// Implementation of the ITS clusterer V2 class //
17// //
18// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
19// //
20///////////////////////////////////////////////////////////////////////////
21
04366a57 22
23#include "AliITSClusterFinderV2SSD.h"
00a7cc50 24#include "AliITSRecPoint.h"
1f3e997f 25#include "AliITSgeomTGeo.h"
7d62fb64 26#include "AliITSDetTypeRec.h"
04366a57 27#include "AliRawReader.h"
28#include "AliITSRawStreamSSD.h"
04366a57 29#include <TClonesArray.h>
04366a57 30#include "AliITSdigitSSD.h"
3a4139a2 31#include "AliITSCalibrationSSD.h"
04366a57 32
c157c94e 33Short_t *AliITSClusterFinderV2SSD::fPairs = 0x0;
34Int_t AliITSClusterFinderV2SSD::fPairsSize = 0;
35
04366a57 36ClassImp(AliITSClusterFinderV2SSD)
37
38
e56160b8 39AliITSClusterFinderV2SSD::AliITSClusterFinderV2SSD(AliITSDetTypeRec* dettyp):AliITSClusterFinderV2(dettyp),
1f3e997f 40fLastSSD1(AliITSgeomTGeo::GetModuleIndex(6,1,1)-1),
e56160b8 41fYpitchSSD(0.0095),
42fHwSSD(3.65),
43fHlSSD(2.00),
44fTanP(0.0275),
45fTanN(0.0075){
04366a57 46
47 //Default constructor
48
04366a57 49}
50
51
52void AliITSClusterFinderV2SSD::FindRawClusters(Int_t mod){
53
54 //Find clusters V2
55 SetModule(mod);
56 FindClustersSSD(fDigits);
57
58}
59
60void AliITSClusterFinderV2SSD::FindClustersSSD(TClonesArray *alldigits) {
61 //------------------------------------------------------------
62 // Actual SSD cluster finder
63 //------------------------------------------------------------
3a4139a2 64 AliITSCalibrationSSD* cal = (AliITSCalibrationSSD*)GetResp(fModule);
65 Float_t gain=0;
66
04366a57 67 Int_t smaxall=alldigits->GetEntriesFast();
68 if (smaxall==0) return;
f7b30404 69 // TObjArray *digits = new TObjArray;
70 TObjArray digits;
04366a57 71 for (Int_t i=0;i<smaxall; i++){
72 AliITSdigitSSD *d=(AliITSdigitSSD*)alldigits->UncheckedAt(i);
3a4139a2 73
74 if(d->IsSideP()) gain = cal->GetGainP(d->GetStripNumber());
75 else gain = cal->GetGainN(d->GetStripNumber());
76
77 Float_t q=gain*d->GetSignal(); // calibration brings mip peaks around 120 (in ADC units)
78 q=cal->ADCToKeV(q); // converts the charge in KeV from ADC units
79 //Float_t q=d->GetSignal()/4.29;// temp. fix (for PID purposed - normalis. to be checked)
6490e1c5 80 d->SetSignal(Int_t(q));
3a4139a2 81
04366a57 82 if (d->GetSignal()<3) continue;
f7b30404 83 digits.AddLast(d);
04366a57 84 }
f7b30404 85 Int_t smax = digits.GetEntriesFast();
04366a57 86 if (smax==0) return;
87
88 const Int_t kMax=1000;
89 Int_t np=0, nn=0;
90 Ali1Dcluster pos[kMax], neg[kMax];
91 Float_t y=0., q=0., qmax=0.;
92 Int_t lab[4]={-2,-2,-2,-2};
93
f7b30404 94 AliITSdigitSSD *d=(AliITSdigitSSD*)digits.UncheckedAt(0);
04366a57 95 q += d->GetSignal();
96 y += d->GetCoord2()*d->GetSignal();
97 qmax=d->GetSignal();
98 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
99 Int_t curr=d->GetCoord2();
100 Int_t flag=d->GetCoord1();
101 Int_t *n=&nn;
102 Ali1Dcluster *c=neg;
103 Int_t nd=1;
104 Int_t milab[10];
105 for (Int_t ilab=0;ilab<10;ilab++){
106 milab[ilab]=-2;
107 }
108 milab[0]=d->GetTrack(0); milab[1]=d->GetTrack(1); milab[2]=d->GetTrack(2);
109
110 for (Int_t s=1; s<smax; s++) {
f7b30404 111 d=(AliITSdigitSSD*)digits.UncheckedAt(s);
04366a57 112 Int_t strip=d->GetCoord2();
113 if ((strip-curr) > 1 || flag!=d->GetCoord1()) {
114 c[*n].SetY(y/q);
115 c[*n].SetQ(q);
116 c[*n].SetNd(nd);
117 CheckLabels2(milab);
118 c[*n].SetLabels(milab);
119 //Split suspiciously big cluster
120 /*
121 if (nd>10&&nd<16){
122 c[*n].SetY(y/q-0.3*nd);
123 c[*n].SetQ(0.5*q);
124 (*n)++;
125 if (*n==MAX) {
126 Error("FindClustersSSD","Too many 1D clusters !");
127 return;
128 }
129 c[*n].SetY(y/q-0.0*nd);
130 c[*n].SetQ(0.5*q);
131 c[*n].SetNd(nd);
132 (*n)++;
133 if (*n==MAX) {
134 Error("FindClustersSSD","Too many 1D clusters !");
135 return;
136 }
137 //
138 c[*n].SetY(y/q+0.3*nd);
139 c[*n].SetQ(0.5*q);
140 c[*n].SetNd(nd);
141 c[*n].SetLabels(milab);
142 }
143 else{
144 */
145 if (nd>4&&nd<25) {
146 c[*n].SetY(y/q-0.25*nd);
147 c[*n].SetQ(0.5*q);
148 (*n)++;
149 if (*n==kMax) {
150 Error("FindClustersSSD","Too many 1D clusters !");
151 return;
152 }
153 c[*n].SetY(y/q+0.25*nd);
154 c[*n].SetQ(0.5*q);
155 c[*n].SetNd(nd);
156 c[*n].SetLabels(milab);
157 }
158 (*n)++;
159 if (*n==kMax) {
160 Error("FindClustersSSD","Too many 1D clusters !");
161 return;
162 }
163 y=q=qmax=0.;
164 nd=0;
165 lab[0]=lab[1]=lab[2]=-2;
166 //
167 for (Int_t ilab=0;ilab<10;ilab++){
168 milab[ilab]=-2;
169 }
170 //
171 if (flag!=d->GetCoord1()) { n=&np; c=pos; }
172 }
173 flag=d->GetCoord1();
174 q += d->GetSignal();
175 y += d->GetCoord2()*d->GetSignal();
176 nd++;
177 if (d->GetSignal()>qmax) {
178 qmax=d->GetSignal();
179 lab[0]=d->GetTrack(0); lab[1]=d->GetTrack(1); lab[2]=d->GetTrack(2);
180 }
181 for (Int_t ilab=0;ilab<10;ilab++) {
182 if (d->GetTrack(ilab)>=0) AddLabel(milab, (d->GetTrack(ilab)));
183 }
184 curr=strip;
185 }
186 c[*n].SetY(y/q);
187 c[*n].SetQ(q);
188 c[*n].SetNd(nd);
189 c[*n].SetLabels(lab);
190 //Split suspiciously big cluster
191 if (nd>4 && nd<25) {
192 c[*n].SetY(y/q-0.25*nd);
193 c[*n].SetQ(0.5*q);
194 (*n)++;
195 if (*n==kMax) {
196 Error("FindClustersSSD","Too many 1D clusters !");
197 return;
198 }
199 c[*n].SetY(y/q+0.25*nd);
200 c[*n].SetQ(0.5*q);
201 c[*n].SetNd(nd);
202 c[*n].SetLabels(lab);
203 }
204 (*n)++;
205 if (*n==kMax) {
206 Error("FindClustersSSD","Too many 1D clusters !");
207 return;
208 }
209
210 FindClustersSSD(neg, nn, pos, np);
211}
212
213
214void AliITSClusterFinderV2SSD::RawdataToClusters(AliRawReader* rawReader,TClonesArray** clusters){
215
216 //------------------------------------------------------------
217 // This function creates ITS clusters from raw data
218 //------------------------------------------------------------
219 rawReader->Reset();
220 AliITSRawStreamSSD inputSSD(rawReader);
221 FindClustersSSD(&inputSSD,clusters);
222
223}
224
3a4139a2 225void AliITSClusterFinderV2SSD::FindClustersSSD(AliITSRawStreamSSD* input,
04366a57 226 TClonesArray** clusters)
227{
228 //------------------------------------------------------------
229 // Actual SSD cluster finder for raw data
230 //------------------------------------------------------------
231 Int_t nClustersSSD = 0;
232 const Int_t kMax = 1000;
233 Ali1Dcluster clusters1D[2][kMax];
234 Int_t nClusters[2] = {0, 0};
235 Int_t lab[3]={-2,-2,-2};
236 Float_t q = 0.;
237 Float_t y = 0.;
238 Int_t nDigits = 0;
239 Int_t prevStrip = -1;
240 Int_t prevFlag = -1;
241 Int_t prevModule = -1;
3a4139a2 242 Float_t gain=0;
b4704be3 243 AliITSCalibrationSSD* cal=NULL;
3a4139a2 244
04366a57 245
246 // read raw data input stream
247 while (kTRUE) {
248 Bool_t next = input->Next();
249
3a4139a2 250 if(input->GetSignal()<(3*4.) && next) continue;
04366a57 251 // check if a new cluster starts
252 Int_t strip = input->GetCoord2();
253 Int_t flag = input->GetCoord1();
254 if ((!next || (input->GetModuleID() != prevModule)||
255 (strip-prevStrip > 1) || (flag != prevFlag)) &&
256 (nDigits > 0)) {
257 if (nClusters[prevFlag] == kMax) {
258 Error("FindClustersSSD", "Too many 1D clusters !");
259 return;
260 }
261 Ali1Dcluster& cluster = clusters1D[prevFlag][nClusters[prevFlag]++];
262 cluster.SetY(y/q);
263 cluster.SetQ(q);
264 cluster.SetNd(nDigits);
265 cluster.SetLabels(lab);
266
267 //Split suspiciously big cluster
268 if (nDigits > 4&&nDigits < 25) {
269 cluster.SetY(y/q - 0.25*nDigits);
270 cluster.SetQ(0.5*q);
271 if (nClusters[prevFlag] == kMax) {
272 Error("FindClustersSSD", "Too many 1D clusters !");
273 return;
274 }
275 Ali1Dcluster& cluster2 = clusters1D[prevFlag][nClusters[prevFlag]++];
276 cluster2.SetY(y/q + 0.25*nDigits);
277 cluster2.SetQ(0.5*q);
278 cluster2.SetNd(nDigits);
279 cluster2.SetLabels(lab);
280 }
281 y = q = 0.;
282 nDigits = 0;
283 }
284
285 if (!next || (input->GetModuleID() != prevModule)) {
286 Int_t iModule = prevModule;
287
288 // when all data from a module was read, search for clusters
289 if (prevFlag >= 0) {
00a7cc50 290 clusters[iModule] = new TClonesArray("AliITSRecPoint");
04366a57 291 fModule = iModule;
292 FindClustersSSD(&clusters1D[0][0], nClusters[0],
293 &clusters1D[1][0], nClusters[1], clusters[iModule]);
294 Int_t nClusters = clusters[iModule]->GetEntriesFast();
295 nClustersSSD += nClusters;
296 }
297
298 if (!next) break;
299 nClusters[0] = nClusters[1] = 0;
300 y = q = 0.;
301 nDigits = 0;
3a4139a2 302
303 cal = (AliITSCalibrationSSD*)GetResp(input->GetModuleID());
304
04366a57 305 }
306
3a4139a2 307 if(input->GetSideFlag()==0) gain = cal->GetGainP(input->GetStrip());
308 else gain = cal->GetGainN(input->GetStrip());
309
04366a57 310 // add digit to current cluster
3a4139a2 311 q += cal->ADCToKeV( gain * input->GetSignal() ); // signal is corrected for gain and converted in KeV
312 y += strip * cal->ADCToKeV( gain * input->GetSignal() );
04366a57 313 nDigits++;
314 prevStrip = strip;
315 prevFlag = flag;
316 prevModule = input->GetModuleID();
317
318 }
319
320 Info("FindClustersSSD", "found clusters in ITS SSD: %d", nClustersSSD);
321}
322
323void AliITSClusterFinderV2SSD::
324FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
325 Ali1Dcluster* pos, Int_t np,
326 TClonesArray *clusters) {
327 //------------------------------------------------------------
328 // Actual SSD cluster finder
329 //------------------------------------------------------------
b4704be3 330
331 const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(fModule);
332
04366a57 333 TClonesArray &cl=*clusters;
334 //
335 Float_t tanp=fTanP, tann=fTanN;
336 if (fModule>fLastSSD1) {tann=fTanP; tanp=fTanN;}
337 Int_t idet=fNdet[fModule];
338 Int_t ncl=0;
339 //
340 Int_t negativepair[30000];
341 Int_t cnegative[3000];
342 Int_t cused1[3000];
343 Int_t positivepair[30000];
344 Int_t cpositive[3000];
345 Int_t cused2[3000];
346 for (Int_t i=0;i<3000;i++) {cnegative[i]=0; cused1[i]=0;}
347 for (Int_t i=0;i<3000;i++) {cpositive[i]=0; cused2[i]=0;}
78538bfe 348 for (Int_t i=0;i<30000;i++) {negativepair[i]=0; positivepair[i]=0;}
c157c94e 349
350 if ((np*nn) > fPairsSize) {
351 if (fPairs) delete [] fPairs;
352 fPairsSize = 4*np*nn;
353 fPairs = new Short_t[fPairsSize];
354 }
355 memset(fPairs,0,sizeof(Short_t)*np*nn);
04366a57 356 //
357 // find available pairs
358 //
359 for (Int_t i=0; i<np; i++) {
360 Float_t yp=pos[i].GetY()*fYpitchSSD;
361 if (pos[i].GetQ()<3) continue;
362 for (Int_t j=0; j<nn; j++) {
363 if (neg[j].GetQ()<3) continue;
364 Float_t yn=neg[j].GetY()*fYpitchSSD;
365 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
366 Float_t yt=yn + tann*zt;
367 zt-=fHlSSD; yt-=fHwSSD;
368 if (TMath::Abs(yt)<fHwSSD+0.01)
369 if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
370 negativepair[i*10+cnegative[i]] =j; //index
371 positivepair[j*10+cpositive[j]] =i;
372 cnegative[i]++; //counters
373 cpositive[j]++;
c157c94e 374 fPairs[i*nn+j]=100;
04366a57 375 }
376 }
377 }
378 //
379 for (Int_t i=0; i<np; i++) {
380 Float_t yp=pos[i].GetY()*fYpitchSSD;
381 if (pos[i].GetQ()<3) continue;
382 for (Int_t j=0; j<nn; j++) {
383 if (neg[j].GetQ()<3) continue;
384 if (cpositive[j]&&cnegative[i]) continue;
385 Float_t yn=neg[j].GetY()*fYpitchSSD;
386 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
387 Float_t yt=yn + tann*zt;
388 zt-=fHlSSD; yt-=fHwSSD;
389 if (TMath::Abs(yt)<fHwSSD+0.1)
390 if (TMath::Abs(zt)<fHlSSD+0.15) {
391 if (cnegative[i]==0) pos[i].SetNd(100); // not available pair
392 if (cpositive[j]==0) neg[j].SetNd(100); // not available pair
393 negativepair[i*10+cnegative[i]] =j; //index
394 positivepair[j*10+cpositive[j]] =i;
395 cnegative[i]++; //counters
396 cpositive[j]++;
c157c94e 397 fPairs[i*nn+j]=100;
04366a57 398 }
399 }
400 }
401 //
402 Float_t lp[5];
403 Int_t milab[10];
404 Double_t ratio;
405
406 //
407 // sign gold tracks
408 //
409 for (Int_t ip=0;ip<np;ip++){
410 Float_t ybest=1000,zbest=1000,qbest=0;
411 //
412 // select gold clusters
413 if ( (cnegative[ip]==1) && cpositive[negativepair[10*ip]]==1){
414 Float_t yp=pos[ip].GetY()*fYpitchSSD;
415 Int_t j = negativepair[10*ip];
416 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
417 //
418 Float_t yn=neg[j].GetY()*fYpitchSSD;
419 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
420 Float_t yt=yn + tann*zt;
421 zt-=fHlSSD; yt-=fHwSSD;
422 ybest=yt; zbest=zt;
423 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
b4704be3 424 {
425 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
426 mT2L->MasterToLocal(loc,trk);
427 lp[0]=trk[1];
428 lp[1]=trk[2];
429 }
04366a57 430 lp[2]=0.0025*0.0025; //SigmaY2
431 lp[3]=0.110*0.110; //SigmaZ2
432
433 lp[4]=qbest; //Q
434 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
435 for (Int_t ilab=0;ilab<3;ilab++){
436 milab[ilab] = pos[ip].GetLabel(ilab);
437 milab[ilab+3] = neg[j].GetLabel(ilab);
438 }
439 //
440 CheckLabels2(milab);
441 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
442 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
00a7cc50 443 AliITSRecPoint * cl2;
444 if(clusters){
75fb37cc 445 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 446 cl2->SetChargeRatio(ratio);
447 cl2->SetType(1);
c157c94e 448 fPairs[ip*nn+j]=1;
00a7cc50 449 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
450 cl2->SetType(2);
c157c94e 451 fPairs[ip*nn+j]=2;
00a7cc50 452 }
453 cused1[ip]++;
454 cused2[j]++;
455
456 }
04366a57 457 else{
75fb37cc 458 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 459 cl2->SetChargeRatio(ratio);
460 cl2->SetType(1);
c157c94e 461 fPairs[ip*nn+j]=1;
00a7cc50 462 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
463 cl2->SetType(2);
c157c94e 464 fPairs[ip*nn+j]=2;
00a7cc50 465 }
466 cused1[ip]++;
467 cused2[j]++;
468 fDetTypeRec->AddRecPoint(*cl2);
04366a57 469 }
470 ncl++;
04366a57 471 }
472 }
473
474 for (Int_t ip=0;ip<np;ip++){
475 Float_t ybest=1000,zbest=1000,qbest=0;
476 //
477 //
478 // select "silber" cluster
479 if ( cnegative[ip]==1 && cpositive[negativepair[10*ip]]==2){
480 Int_t in = negativepair[10*ip];
481 Int_t ip2 = positivepair[10*in];
482 if (ip2==ip) ip2 = positivepair[10*in+1];
483 Float_t pcharge = pos[ip].GetQ()+pos[ip2].GetQ();
484 if (TMath::Abs(pcharge-neg[in].GetQ())<10){
485 //
486 // add first pair
c157c94e 487 if (fPairs[ip*nn+in]==100){ //
04366a57 488 Float_t yp=pos[ip].GetY()*fYpitchSSD;
489 Float_t yn=neg[in].GetY()*fYpitchSSD;
490 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
491 Float_t yt=yn + tann*zt;
492 zt-=fHlSSD; yt-=fHwSSD;
493 ybest =yt; zbest=zt;
494 qbest =pos[ip].GetQ();
b4704be3 495 {
496 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
497 mT2L->MasterToLocal(loc,trk);
498 lp[0]=trk[1];
499 lp[1]=trk[2];
500 }
04366a57 501 lp[2]=0.0025*0.0025; //SigmaY2
502 lp[3]=0.110*0.110; //SigmaZ2
503
504 lp[4]=qbest; //Q
505 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
506 for (Int_t ilab=0;ilab<3;ilab++){
507 milab[ilab] = pos[ip].GetLabel(ilab);
508 milab[ilab+3] = neg[in].GetLabel(ilab);
509 }
510 //
511 CheckLabels2(milab);
512 ratio = (pos[ip].GetQ()-neg[in].GetQ())/(pos[ip].GetQ()+neg[in].GetQ());
513 milab[3]=(((ip<<10) + in)<<10) + idet; // pos|neg|det
514 Int_t info[3] = {pos[ip].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
515
00a7cc50 516 AliITSRecPoint * cl2;
517 if(clusters){
75fb37cc 518 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 519 cl2->SetChargeRatio(ratio);
520 cl2->SetType(5);
c157c94e 521 fPairs[ip*nn+in] = 5;
00a7cc50 522 if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
523 cl2->SetType(6);
c157c94e 524 fPairs[ip*nn+in] = 6;
00a7cc50 525 }
526 }
04366a57 527 else{
75fb37cc 528 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 529 cl2->SetChargeRatio(ratio);
530 cl2->SetType(5);
c157c94e 531 fPairs[ip*nn+in] = 5;
00a7cc50 532 if ((pos[ip].GetNd()+neg[in].GetNd())>6){ //multi cluster
533 cl2->SetType(6);
c157c94e 534 fPairs[ip*nn+in] = 6;
00a7cc50 535 }
536
537 fDetTypeRec->AddRecPoint(*cl2);
04366a57 538 }
539 ncl++;
04366a57 540 }
00a7cc50 541
04366a57 542 //
543 // add second pair
544
00a7cc50 545 // if (!(cused1[ip2] || cused2[in])){ //
c157c94e 546 if (fPairs[ip2*nn+in]==100){
04366a57 547 Float_t yp=pos[ip2].GetY()*fYpitchSSD;
548 Float_t yn=neg[in].GetY()*fYpitchSSD;
549 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
550 Float_t yt=yn + tann*zt;
551 zt-=fHlSSD; yt-=fHwSSD;
552 ybest =yt; zbest=zt;
553 qbest =pos[ip2].GetQ();
b4704be3 554 {
555 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
556 mT2L->MasterToLocal(loc,trk);
557 lp[0]=trk[1];
558 lp[1]=trk[2];
559 }
04366a57 560 lp[2]=0.0025*0.0025; //SigmaY2
561 lp[3]=0.110*0.110; //SigmaZ2
562
563 lp[4]=qbest; //Q
564 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
565 for (Int_t ilab=0;ilab<3;ilab++){
566 milab[ilab] = pos[ip2].GetLabel(ilab);
567 milab[ilab+3] = neg[in].GetLabel(ilab);
568 }
569 //
570 CheckLabels2(milab);
571 ratio = (pos[ip2].GetQ()-neg[in].GetQ())/(pos[ip2].GetQ()+neg[in].GetQ());
572 milab[3]=(((ip2<<10) + in)<<10) + idet; // pos|neg|det
573 Int_t info[3] = {pos[ip2].GetNd(),neg[in].GetNd(),fNlayer[fModule]};
00a7cc50 574
575 AliITSRecPoint * cl2;
576 if(clusters){
75fb37cc 577 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 578 cl2->SetChargeRatio(ratio);
579 cl2->SetType(5);
c157c94e 580 fPairs[ip2*nn+in] =5;
00a7cc50 581 if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
582 cl2->SetType(6);
c157c94e 583 fPairs[ip2*nn+in] =6;
00a7cc50 584 }
585 }
04366a57 586 else{
75fb37cc 587 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 588 cl2->SetChargeRatio(ratio);
589 cl2->SetType(5);
c157c94e 590 fPairs[ip2*nn+in] =5;
00a7cc50 591 if ((pos[ip2].GetNd()+neg[in].GetNd())>6){ //multi cluster
592 cl2->SetType(6);
c157c94e 593 fPairs[ip2*nn+in] =6;
00a7cc50 594 }
595
596 fDetTypeRec->AddRecPoint(*cl2);
04366a57 597 }
598 ncl++;
04366a57 599 }
600 cused1[ip]++;
601 cused1[ip2]++;
602 cused2[in]++;
603 }
604 }
605 }
00a7cc50 606
04366a57 607 //
608 for (Int_t jn=0;jn<nn;jn++){
609 if (cused2[jn]) continue;
610 Float_t ybest=1000,zbest=1000,qbest=0;
611 // select "silber" cluster
612 if ( cpositive[jn]==1 && cnegative[positivepair[10*jn]]==2){
613 Int_t ip = positivepair[10*jn];
614 Int_t jn2 = negativepair[10*ip];
615 if (jn2==jn) jn2 = negativepair[10*ip+1];
616 Float_t pcharge = neg[jn].GetQ()+neg[jn2].GetQ();
617 //
618 if (TMath::Abs(pcharge-pos[ip].GetQ())<10){
619 //
620 // add first pair
621 // if (!(cused1[ip]||cused2[jn])){
c157c94e 622 if (fPairs[ip*nn+jn]==100){
04366a57 623 Float_t yn=neg[jn].GetY()*fYpitchSSD;
624 Float_t yp=pos[ip].GetY()*fYpitchSSD;
625 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
626 Float_t yt=yn + tann*zt;
627 zt-=fHlSSD; yt-=fHwSSD;
628 ybest =yt; zbest=zt;
629 qbest =neg[jn].GetQ();
b4704be3 630 {
631 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
632 mT2L->MasterToLocal(loc,trk);
633 lp[0]=trk[1];
634 lp[1]=trk[2];
635 }
04366a57 636 lp[2]=0.0025*0.0025; //SigmaY2
637 lp[3]=0.110*0.110; //SigmaZ2
638
639 lp[4]=qbest; //Q
640 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
641 for (Int_t ilab=0;ilab<3;ilab++){
642 milab[ilab] = pos[ip].GetLabel(ilab);
643 milab[ilab+3] = neg[jn].GetLabel(ilab);
644 }
645 //
646 CheckLabels2(milab);
647 ratio = (pos[ip].GetQ()-neg[jn].GetQ())/(pos[ip].GetQ()+neg[jn].GetQ());
648 milab[3]=(((ip<<10) + jn)<<10) + idet; // pos|neg|det
649 Int_t info[3] = {pos[ip].GetNd(),neg[jn].GetNd(),fNlayer[fModule]};
650
00a7cc50 651 AliITSRecPoint * cl2;
652 if(clusters){
75fb37cc 653 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 654 cl2->SetChargeRatio(ratio);
655 cl2->SetType(7);
c157c94e 656 fPairs[ip*nn+jn] =7;
00a7cc50 657 if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
658 cl2->SetType(8);
c157c94e 659 fPairs[ip*nn+jn]=8;
00a7cc50 660 }
661
662 }
04366a57 663 else{
75fb37cc 664 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 665 cl2->SetChargeRatio(ratio);
666 cl2->SetType(7);
c157c94e 667 fPairs[ip*nn+jn] =7;
00a7cc50 668 if ((pos[ip].GetNd()+neg[jn].GetNd())>6){ //multi cluster
669 cl2->SetType(8);
c157c94e 670 fPairs[ip*nn+jn]=8;
00a7cc50 671 }
672
673 fDetTypeRec->AddRecPoint(*cl2);
04366a57 674 }
675 ncl++;
04366a57 676 }
677 //
678 // add second pair
679 // if (!(cused1[ip]||cused2[jn2])){
c157c94e 680 if (fPairs[ip*nn+jn2]==100){
04366a57 681 Float_t yn=neg[jn2].GetY()*fYpitchSSD;
682 Double_t yp=pos[ip].GetY()*fYpitchSSD;
683 Double_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
684 Double_t yt=yn + tann*zt;
685 zt-=fHlSSD; yt-=fHwSSD;
686 ybest =yt; zbest=zt;
687 qbest =neg[jn2].GetQ();
b4704be3 688 {
689 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
690 mT2L->MasterToLocal(loc,trk);
691 lp[0]=trk[1];
692 lp[1]=trk[2];
693 }
04366a57 694 lp[2]=0.0025*0.0025; //SigmaY2
695 lp[3]=0.110*0.110; //SigmaZ2
696
697 lp[4]=qbest; //Q
698 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
699 for (Int_t ilab=0;ilab<3;ilab++){
700 milab[ilab] = pos[ip].GetLabel(ilab);
701 milab[ilab+3] = neg[jn2].GetLabel(ilab);
702 }
703 //
704 CheckLabels2(milab);
705 ratio = (pos[ip].GetQ()-neg[jn2].GetQ())/(pos[ip].GetQ()+neg[jn2].GetQ());
706 milab[3]=(((ip<<10) + jn2)<<10) + idet; // pos|neg|det
707 Int_t info[3] = {pos[ip].GetNd(),neg[jn2].GetNd(),fNlayer[fModule]};
00a7cc50 708 AliITSRecPoint * cl2;
709 if(clusters){
75fb37cc 710 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 711 cl2->SetChargeRatio(ratio);
c157c94e 712 fPairs[ip*nn+jn2]=7;
00a7cc50 713 cl2->SetType(7);
714 if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
715 cl2->SetType(8);
c157c94e 716 fPairs[ip*nn+jn2]=8;
00a7cc50 717 }
718
719 }
04366a57 720 else{
75fb37cc 721 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 722 cl2->SetChargeRatio(ratio);
c157c94e 723 fPairs[ip*nn+jn2]=7;
00a7cc50 724 cl2->SetType(7);
725 if ((pos[ip].GetNd()+neg[jn2].GetNd())>6){ //multi cluster
726 cl2->SetType(8);
c157c94e 727 fPairs[ip*nn+jn2]=8;
00a7cc50 728 }
729
730 fDetTypeRec->AddRecPoint(*cl2);
04366a57 731 }
732
733 ncl++;
04366a57 734 }
735 cused1[ip]++;
736 cused2[jn]++;
737 cused2[jn2]++;
738 }
739 }
740 }
741
742 for (Int_t ip=0;ip<np;ip++){
743 Float_t ybest=1000,zbest=1000,qbest=0;
744 //
745 // 2x2 clusters
746 //
747 if ( (cnegative[ip]<5) && cpositive[negativepair[10*ip]]<5){
748 Float_t minchargediff =4.;
749 Int_t j=-1;
750 for (Int_t di=0;di<cnegative[ip];di++){
751 Int_t jc = negativepair[ip*10+di];
752 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
753 if (TMath::Abs(chargedif)<minchargediff){
754 j =jc;
755 minchargediff = TMath::Abs(chargedif);
756 }
757 }
758 if (j<0) continue; // not proper cluster
759 Int_t count =0;
760 for (Int_t di=0;di<cnegative[ip];di++){
761 Int_t jc = negativepair[ip*10+di];
762 Float_t chargedif = pos[ip].GetQ()-neg[jc].GetQ();
763 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
764 }
765 if (count>1) continue; // more than one "proper" cluster for positive
766 //
767 count =0;
768 for (Int_t dj=0;dj<cpositive[j];dj++){
769 Int_t ic = positivepair[j*10+dj];
770 Float_t chargedif = pos[ic].GetQ()-neg[j].GetQ();
771 if (TMath::Abs(chargedif)<minchargediff+3.) count++;
772 }
773 if (count>1) continue; // more than one "proper" cluster for negative
774
775 Int_t jp = 0;
776
777 count =0;
778 for (Int_t dj=0;dj<cnegative[jp];dj++){
779 Int_t ic = positivepair[jp*10+dj];
780 Float_t chargedif = pos[ic].GetQ()-neg[jp].GetQ();
781 if (TMath::Abs(chargedif)<minchargediff+4.) count++;
782 }
783 if (count>1) continue;
c157c94e 784 if (fPairs[ip*nn+j]<100) continue;
04366a57 785 //
786 //almost gold clusters
787 Float_t yp=pos[ip].GetY()*fYpitchSSD;
788 Float_t yn=neg[j].GetY()*fYpitchSSD;
789 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
790 Float_t yt=yn + tann*zt;
791 zt-=fHlSSD; yt-=fHwSSD;
792 ybest=yt; zbest=zt;
793 qbest=0.5*(pos[ip].GetQ()+neg[j].GetQ());
b4704be3 794 {
795 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
796 mT2L->MasterToLocal(loc,trk);
797 lp[0]=trk[1];
798 lp[1]=trk[2];
799 }
04366a57 800 lp[2]=0.0025*0.0025; //SigmaY2
801 lp[3]=0.110*0.110; //SigmaZ2
802 lp[4]=qbest; //Q
803 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
804 for (Int_t ilab=0;ilab<3;ilab++){
805 milab[ilab] = pos[ip].GetLabel(ilab);
806 milab[ilab+3] = neg[j].GetLabel(ilab);
807 }
808 //
809 CheckLabels2(milab);
810 ratio = (pos[ip].GetQ()-neg[j].GetQ())/(pos[ip].GetQ()+neg[j].GetQ());
811 milab[3]=(((ip<<10) + j)<<10) + idet; // pos|neg|det
812 Int_t info[3] = {pos[ip].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
00a7cc50 813 AliITSRecPoint * cl2;
814 if(clusters){
75fb37cc 815 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 816 cl2->SetChargeRatio(ratio);
817 cl2->SetType(10);
c157c94e 818 fPairs[ip*nn+j]=10;
00a7cc50 819 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
820 cl2->SetType(11);
c157c94e 821 fPairs[ip*nn+j]=11;
00a7cc50 822 }
823 cused1[ip]++;
824 cused2[j]++;
04366a57 825 }
00a7cc50 826 else{
75fb37cc 827 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 828 cl2->SetChargeRatio(ratio);
829 cl2->SetType(10);
c157c94e 830 fPairs[ip*nn+j]=10;
00a7cc50 831 if ((pos[ip].GetNd()+neg[j].GetNd())>6){ //multi cluster
832 cl2->SetType(11);
c157c94e 833 fPairs[ip*nn+j]=11;
00a7cc50 834 }
835 cused1[ip]++;
836 cused2[j]++;
837
838 fDetTypeRec->AddRecPoint(*cl2);
839 }
04366a57 840 ncl++;
04366a57 841 }
842
843 }
844
845 //
846 for (Int_t i=0; i<np; i++) {
847 Float_t ybest=1000,zbest=1000,qbest=0;
848 Float_t yp=pos[i].GetY()*fYpitchSSD;
849 if (pos[i].GetQ()<3) continue;
850 for (Int_t j=0; j<nn; j++) {
851 // for (Int_t di = 0;di<cpositive[i];di++){
852 // Int_t j = negativepair[10*i+di];
853 if (neg[j].GetQ()<3) continue;
854 if (cused2[j]||cused1[i]) continue;
c157c94e 855 if (fPairs[i*nn+j]>0 &&fPairs[i*nn+j]<100) continue;
04366a57 856 ratio = (pos[i].GetQ()-neg[j].GetQ())/(pos[i].GetQ()+neg[j].GetQ());
857 Float_t yn=neg[j].GetY()*fYpitchSSD;
858 Float_t zt=(2*fHlSSD*tanp + yp - yn)/(tann+tanp);
859 Float_t yt=yn + tann*zt;
860 zt-=fHlSSD; yt-=fHwSSD;
861 if (TMath::Abs(yt)<fHwSSD+0.01)
862 if (TMath::Abs(zt)<fHlSSD+0.01*(neg[j].GetNd()+pos[i].GetNd())) {
863 ybest=yt; zbest=zt;
864 qbest=0.5*(pos[i].GetQ()+neg[j].GetQ());
b4704be3 865 {
866 Double_t loc[3]={ybest,0.,zbest},trk[3]={0.,0.,0.};
867 mT2L->MasterToLocal(loc,trk);
868 lp[0]=trk[1];
869 lp[1]=trk[2];
870 }
04366a57 871 lp[2]=0.0025*0.0025; //SigmaY2
872 lp[3]=0.110*0.110; //SigmaZ2
873
874 lp[4]=qbest; //Q
875 for (Int_t ilab=0;ilab<10;ilab++) milab[ilab]=-2;
876 for (Int_t ilab=0;ilab<3;ilab++){
877 milab[ilab] = pos[i].GetLabel(ilab);
878 milab[ilab+3] = neg[j].GetLabel(ilab);
879 }
880 //
881 CheckLabels2(milab);
882 milab[3]=(((i<<10) + j)<<10) + idet; // pos|neg|det
883 Int_t info[3] = {pos[i].GetNd(),neg[j].GetNd(),fNlayer[fModule]};
00a7cc50 884 AliITSRecPoint * cl2;
885 if(clusters){
75fb37cc 886 cl2 = new (cl[ncl]) AliITSRecPoint(milab,lp,info);
00a7cc50 887 cl2->SetChargeRatio(ratio);
888 cl2->SetType(100+cpositive[j]+cnegative[i]);
889 }
04366a57 890 else{
75fb37cc 891 cl2 = new AliITSRecPoint(milab,lp,info);
00a7cc50 892 cl2->SetChargeRatio(ratio);
893 cl2->SetType(100+cpositive[j]+cnegative[i]);
894 fDetTypeRec->AddRecPoint(*cl2);
04366a57 895 }
896 ncl++;
04366a57 897 //cl2->SetType(0);
898 /*
c157c94e 899 if (fPairs[i*nn+j]<100){
900 printf("problem:- %d\n", fPairs[i*nn+j]);
04366a57 901 }
902 if (cnegative[i]<2&&cpositive[j]<2){
c157c94e 903 printf("problem:- %d\n", fPairs[i*nn+j]);
04366a57 904 }
905 */
906 }
907 }
908 }
909
04366a57 910}
911
912
913