Release version of ITS code
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderSSD.cxx
CommitLineData
b0f5e3fc 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/*
17
18Adding rekonstruction facilities
19Piotr Krzysztof Skowronski
20December 1999.
21*/
22
23/*
2415 -18 V 2000
25Eroor counting routines
26Automatic combination routines improved (traps)
27
28*/
29
30#include "AliRun.h"
e8189707 31#include "AliITS.h"
32#include "AliITSMapA1.h"
b0f5e3fc 33#include "AliITSClusterFinderSSD.h"
e8189707 34#include "AliITSclusterSSD.h"
35#include "AliITSpackageSSD.h"
b0f5e3fc 36
37
38const Int_t debug=0;
39
40ClassImp(AliITSClusterFinderSSD)
41
42//____________________________________________________________________
43//
44// Constructor
45//____________________________________________________________________
46//
47
48
49AliITSClusterFinderSSD::AliITSClusterFinderSSD(AliITSsegmentation *seg, TClonesArray *digits, TClonesArray *recp)
50{
51
52 fSegmentation=seg;
53 fDigits=digits;
54 fRecPoints=recp;
55
e8189707 56 fMap = new AliITSMapA1(fSegmentation,fDigits);
57
b0f5e3fc 58 fITS=(AliITS*)gAlice->GetModule("ITS");
59
60 fClusterP = new TClonesArray ("AliITSclusterSSD",200);
61 fNClusterP =0;
62
63 fClusterN= new TClonesArray ("AliITSclusterSSD",200);
64 fNClusterN =0;
65
66 fPackages = new TClonesArray ("AliITSpackageSSD",200); //packages
67 fNPackages = 0;
68
69
70 fDigitsIndexP = new TArrayI(300);
71 fNDigitsP = 0;
72
73 fDigitsIndexN = new TArrayI(300);
74 fNDigitsN = 0;
75
76 SetAlpha1(1000);
77 SetAlpha2(1000);
78 SetAlpha3(1000);
79
80
81 fPitch = fSegmentation->Dpx(0);
82 Float_t StereoP,StereoN;
83 fSegmentation->Angles(StereoP,StereoN);
84 fTanP=TMath::Tan(StereoP);
85 fTanN=TMath::Tan(StereoN);
e8189707 86
be33dccb 87 fPNsignalRatio=7./8.; // warning: hard-wired number
e8189707 88
b0f5e3fc 89}
90
91//-------------------------------------------------------
92AliITSClusterFinderSSD::~AliITSClusterFinderSSD() {
93
94 delete fClusterP;
95 delete fClusterN;
96 delete fPackages;
97 delete fDigitsIndexP;
98 delete fDigitsIndexN;
e8189707 99
100 delete fMap;
101
b0f5e3fc 102}
103
104//-------------------------------------------------------
105void AliITSClusterFinderSSD::InitReconstruction()
106{
107
108 register Int_t i; //iterator
109
e8189707 110 for (i=0;i<fNClusterP;i++)
b0f5e3fc 111 {
112 fClusterP->RemoveAt(i);
113 }
114 fNClusterP =0;
e8189707 115 for (i=0;i<fNClusterN;i++)
b0f5e3fc 116 {
117 fClusterN->RemoveAt(i);
118 }
119 fNClusterN=0;
120
e8189707 121 for (i=0;i<fNPackages;i++)
b0f5e3fc 122 {
123 fPackages->RemoveAt(i);
124 }
125
126 fNPackages = 0;
127 fNDigitsP=0;
128 fNDigitsN=0;
129
130 Float_t StereoP,StereoN;
131 fSegmentation->Angles(StereoP,StereoN);
132
133 CalcStepFactor (StereoP,StereoN);
134
135 if (debug) cout<<"fSFF = "<<fSFF<<" fSFB = "<<fSFB<<"\n";
136}
137
138
139//---------------------------------------------
140void AliITSClusterFinderSSD::FindRawClusters()
141{
b0f5e3fc 142//Piotr Krzysztof Skowronski
143//Warsaw University of Technology
144//skowron@if.pw.edu.pl
145
146// This function findes out all clusters belonging to one module
147// 1. Zeroes all space after previous module reconstruction
148// 2. Finds all neighbouring digits
149// 3. If necesery, resolves for each group of neighbouring digits
150// how many clusters creates it.
151// 4. Creates packages
152// 5. Creates clusters
153
154 InitReconstruction(); //ad. 1
e8189707 155 fMap->FillMap();
b0f5e3fc 156 FillDigitsIndex();
157 SortDigits();
158 FindNeighbouringDigits(); //ad. 2
159 SeparateOverlappedClusters(); //ad. 3
160 ClustersToPackages(); //ad. 4
161 ConsumeClusters();
162 PackagesToPoints(); //ad. 5
163 ReconstructNotConsumedClusters();
e8189707 164
165 fMap->ClearMap();
b0f5e3fc 166}
167
168
169//-------------------------------------------------
170void AliITSClusterFinderSSD::FindNeighbouringDigits()
171{
172
173
174//Piotr Krzysztof Skowronski
175//Warsaw University of Technology
176//skowron@if.pw.edu.pl
177
178 register Int_t i;
179
180 //If there are any digits on this side, create 1st Cluster,
181 // add to it this digit, and increment number of clusters
182
183 if ((fNDigitsP>0 ) && (fNDigitsN > 0 )) {
184
185 Int_t currentstripNo;
186 Int_t *dbuffer = new Int_t [300]; //buffer for strip numbers
187 Int_t dnumber; //curent number of digits in buffer
188 TArrayI &lDigitsIndexP = *fDigitsIndexP;
189 TArrayI &lDigitsIndexN = *fDigitsIndexN;
190 TObjArray &lDigits=*fDigits;
191 TClonesArray &lClusterP = *fClusterP;
192 TClonesArray &lClusterN = *fClusterN;
193
194 //process P side
195 dnumber = 1;
196 dbuffer[0]=lDigitsIndexP[0];
197 //If next digit is a neighbour of previous, adds to last cluster this digit
e8189707 198 for (i=1; i<fNDigitsP; i++) {
b0f5e3fc 199 //reads new digit
200 currentstripNo = ((AliITSdigitSSD*)lDigits[lDigitsIndexP[i]])->
201 GetStripNumber();
202 //check if it is a neighbour of a previous one
203 if ( (((AliITSdigitSSD*)lDigits[lDigitsIndexP[i-1]])->GetStripNumber())
204 == (currentstripNo-1) ) dbuffer[dnumber++]=lDigitsIndexP[i];
205 else {
206 //create a new one side cluster
e8189707 207 new(lClusterP[fNClusterP++]) AliITSclusterSSD(dnumber,dbuffer,fDigits,fgkSIDEP);
b0f5e3fc 208 dbuffer[0]=lDigitsIndexP[i];
209 dnumber = 1;
210 }
211 } // end loop over fNDigitsP
e8189707 212 new(lClusterP[fNClusterP++]) AliITSclusterSSD(dnumber,dbuffer,fDigits,fgkSIDEP);
b0f5e3fc 213
214
215 //process N side
216 //for comments, see above
217 dnumber = 1;
218 dbuffer[0]=lDigitsIndexN[0];
219 //If next digit is a neighbour of previous, adds to last cluster this digit
e8189707 220 for (i=1; i<fNDigitsN; i++) {
b0f5e3fc 221 currentstripNo = ((AliITSdigitSSD*)(lDigits[lDigitsIndexN[i]]))->
222 GetStripNumber();
223 if ( (((AliITSdigitSSD*)lDigits[lDigitsIndexN[i-1]])->GetStripNumber())
224 == (currentstripNo-1) ) dbuffer[dnumber++]=lDigitsIndexN[i];
225 else {
e8189707 226 new(lClusterN[fNClusterN++]) AliITSclusterSSD(dnumber,dbuffer,fDigits,fgkSIDEN);
b0f5e3fc 227 dbuffer[0]=lDigitsIndexN[i];
228 dnumber = 1;
229 }
230 } // end loop over fNDigitsN
e8189707 231 new(lClusterN[fNClusterN++]) AliITSclusterSSD(dnumber,dbuffer,fDigits,fgkSIDEN);
b0f5e3fc 232 delete [] dbuffer;
233
234 } // end condition on NDigits
235
236 if (debug) cout<<"\n Found clusters: fNClusterP = "<<fNClusterP<<" fNClusterN ="<<fNClusterN<<"\n";
237
238}
239/**********************************************************************/
240
241
242void AliITSClusterFinderSSD::SeparateOverlappedClusters()
243{
244//************************************************
245//Piotr Krzysztof Skowronski
246//Warsaw University of Technology
247//skowron@if.pw.edu.pl
248
e8189707 249 register Int_t i; //iterator
b0f5e3fc 250
251 Float_t factor=0.75; // How many percent must be lower signal
252 // on the middle one digit
253 // from its neighbours
254 Int_t signal0; //signal on the strip before the current one
255 Int_t signal1; //signal on the current one signal
256 Int_t signal2; //signal on the strip after the current one
257 TArrayI *splitlist; // List of splits
258 Int_t numerofsplits=0; // number of splits
259 Int_t initPsize = fNClusterP; //initial size of the arrays
260 Int_t initNsize = fNClusterN; //we have to keep it because it will grow
261 // in this function and it doasn't make
262 // sense to pass through it again
263
264 splitlist = new TArrayI(300);
265
e8189707 266 for (i=0;i<initPsize;i++)
b0f5e3fc 267 {
268 if (( ((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits())==1) continue;
269 if (( ((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits())==2) continue;
270 Int_t nj=(((AliITSclusterSSD*)(*fClusterP)[i])->GetNumOfDigits()-1);
e8189707 271 for (Int_t j=1; j<nj; j++)
b0f5e3fc 272 {
273 signal1=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j);
274 signal0=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j-1);
275 signal2=((AliITSclusterSSD*)(*fClusterP)[i])->GetDigitSignal(j+1);
276 //if signal is less then factor*signal of its neighbours
277 if ( (signal1<(factor*signal0)) && (signal1<(factor*signal2)) )
278 {
279 (*splitlist)[numerofsplits++]=j;
280 }
281 } // end loop over number of digits
282 //split this cluster if necessary
e8189707 283 if(numerofsplits>0) SplitCluster(splitlist,numerofsplits,i,fgkSIDEP);
b0f5e3fc 284 numerofsplits=0;
285
286 //in signed places (splitlist)
287 } // end loop over clusters on Pside
288
e8189707 289 for (i=0;i<initNsize;i++) {
b0f5e3fc 290 if (( ((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits())==1) continue;
291 if (( ((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits())==2) continue;
292 Int_t nj=(((AliITSclusterSSD*)(*fClusterN)[i])->GetNumOfDigits()-1);
e8189707 293 for (Int_t j=1; j<nj; j++)
b0f5e3fc 294 {
295 signal1=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j);
296 signal0=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j-1);
297 signal2=((AliITSclusterSSD*)(*fClusterN)[i])->GetDigitSignal(j+1);
298 //if signal is less then factor*signal of its neighbours
299 if ( (signal1<(factor*signal0)) && (signal1<(factor*signal2)) )
300 (*splitlist)[numerofsplits++]=j;
301 } // end loop over number of digits
302 //split this cluster into more clusters
e8189707 303 if(numerofsplits>0) SplitCluster(splitlist,numerofsplits,i,fgkSIDEN);
b0f5e3fc 304 numerofsplits=0; //in signed places (splitlist)
305 } // end loop over clusters on Nside
306
307 delete splitlist;
308}
309
310//-------------------------------------------------------
311void AliITSClusterFinderSSD::SplitCluster(TArrayI *list, Int_t nsplits, Int_t index, Bool_t side)
312{
313
314
315//Piotr Krzysztof Skowronski
316//Warsaw University of Technology
317//skowron@if.pw.edu.pl
318
319 //This function splits one side cluster into more clusters
320 //number of splits is defined by "nsplits"
321 //Place of splits are defined in the TArray "list"
322
323 // For further optimisation: Replace this function by two
324 // specialised ones (each for one side)
325 // save one "if"
326
327 //For comlete comments see AliITSclusterSSD::SplitCluster
328
329
330 register Int_t i; //iterator
331
332 AliITSclusterSSD* curentcluster;
333 Int_t *tmpdigits = new Int_t[100];
334 Int_t NN;
335 // side true means P side
336 if (side) {
337 curentcluster =((AliITSclusterSSD*)((*fClusterP)[index])) ;
e8189707 338 for (i = nsplits; i>0 ;i--) {
b0f5e3fc 339 NN=curentcluster->SplitCluster((*list)[(i-1)],tmpdigits);
340 new ((*fClusterP)[fNClusterP]) AliITSclusterSSD(NN,tmpdigits,fDigits,side);
341 ( (AliITSclusterSSD*)((*fClusterP)[fNClusterP]) )->
342 SetLeftNeighbour(kTRUE);
343 //if left cluster had neighbour on the right before split
344 //new should have it too
345 if ( curentcluster->GetRightNeighbour() )
346 ( (AliITSclusterSSD*)((*fClusterP)[fNClusterP]) )->
347 SetRightNeighbour(kTRUE);
348 else curentcluster->SetRightNeighbour(kTRUE);
349 fNClusterP++;
350 } // end loop over nplits
351 } else {
352 curentcluster =((AliITSclusterSSD*)((*fClusterN)[index]));
e8189707 353 for (i = nsplits; i>0 ;i--) {
b0f5e3fc 354 NN=curentcluster->SplitCluster((*list)[(i-1)],tmpdigits);
355 new ((*fClusterN)[fNClusterN]) AliITSclusterSSD(NN,tmpdigits,fDigits,side);
356 ((AliITSclusterSSD*)((*fClusterN)[fNClusterN]))->
357 SetRightNeighbour(kTRUE);
358 if (curentcluster->GetRightNeighbour())
359 ( (AliITSclusterSSD*)( (*fClusterN)[fNClusterN]) )->
360 SetRightNeighbour(kTRUE);
361 else curentcluster->SetRightNeighbour(kTRUE);
362 fNClusterN++;
363 } // end loop over nplits
364 } // end if side
365 delete []tmpdigits;
366
367}
368
369
370//-------------------------------------------------
371Int_t AliITSClusterFinderSSD::SortDigitsP(Int_t start, Int_t end)
372{
373
374
375//Piotr Krzysztof Skowronski
376//Warsaw University of Technology
377//skowron@if.pw.edu.pl
378
379 Int_t right;
380 Int_t left;
381 if (start != (end - 1) )
382 {
383 left=this->SortDigitsP(start,(start+end)/2);
384 right=this->SortDigitsP((start+end)/2,end);
385 return (left || right);
386 }
387 else
388 {
389 left = ((AliITSdigitSSD*)((*fDigits)[(*fDigitsIndexP)[start]]))->GetStripNumber();
390 right= ((AliITSdigitSSD*)((*fDigits)[(*fDigitsIndexP)[end]]))->GetStripNumber();
391 if( left > right )
392 {
393 Int_t tmp = (*fDigitsIndexP)[start];
394 (*fDigitsIndexP)[start]=(*fDigitsIndexP)[end];
395 (*fDigitsIndexP)[end]=tmp;
396 return 1;
397 }
398 else return 0;
399 }
400}
401
402
403//--------------------------------------------------
404
405Int_t AliITSClusterFinderSSD::SortDigitsN(Int_t start, Int_t end)
406{
407
408
409//Piotr Krzysztof Skowronski
410//Warsaw University of Technology
411//skowron@if.pw.edu.pl
412
413 Int_t right;
414 Int_t left;
415 if (start != (end - 1))
416 {
417 left=this->SortDigitsN(start,(start+end)/2);
418 right=this->SortDigitsN((start+end)/2,end);
419 return (left || right);
420 }
421 else
422 {
423 left =((AliITSdigitSSD*)((*fDigits)[(*fDigitsIndexN)[start]]))->GetStripNumber();
424 right=((AliITSdigitSSD*)((*fDigits)[(*fDigitsIndexN)[end]]))->GetStripNumber();
425 if ( left > right )
426 {
427 Int_t tmp = (*fDigitsIndexN)[start];
428 (*fDigitsIndexN)[start]=(*fDigitsIndexN)[end];
429 (*fDigitsIndexN)[end]=tmp;
430 return 1;
431 }else return 0;
432 }
433}
434
435
436//------------------------------------------------
437void AliITSClusterFinderSSD::FillDigitsIndex()
438{
439
440
441//Piotr Krzysztof Skowronski
442//Warsaw University of Technology
443//skowron@if.pw.edu.pl
444
445 //Fill the indexes of the clusters belonging to a given ITS module
446 //Created by Piotr K. Skowronski, August 7 1999
447
448 Int_t PNs=0, NNs=0;
449 Int_t tmp,bit,k;
450 Int_t N;
451 Int_t i;
452
453 N = fDigits->GetEntriesFast();
454
455 Int_t* PSidx = new Int_t [N*sizeof(Int_t)];
456 Int_t* NSidx = new Int_t [N*sizeof(Int_t)];
457 if (fDigitsIndexP==NULL) fDigitsIndexP = new TArrayI(N);
458 if (fDigitsIndexN==NULL) fDigitsIndexN = new TArrayI(N);
459
460 AliITSdigitSSD *dig;
461
e8189707 462 for ( i = 0 ; i< N; i++ ) {
b0f5e3fc 463 dig=(AliITSdigitSSD*)fDigits->UncheckedAt(i);
464 if(dig->IsSideP()) {
465 bit=1;
466 tmp=dig->GetStripNumber();
467 // I find this totally unnecessary - it's just a
468 // CPU consuming double check
e8189707 469 for( k=0;k<PNs;k++)
b0f5e3fc 470 {
471 if (tmp==PSidx[k])
472 {
473 if (debug) cout<<"Such a digit exists \n";
474 bit=0;
475 }
476 }
477 // end comment
478 if(bit) {
479 fDigitsIndexP->AddAt(i,fNDigitsP++);
480 PSidx[PNs++]=tmp;
481 }
482 } else {
483 bit=1;
484 tmp=dig->GetStripNumber();
485 // same as above
e8189707 486 for( k=0;k<NNs;k++)
b0f5e3fc 487 {
488 if (tmp==NSidx[k])
489 {
490 if (debug) cout<<"Such a digit exists \n";
491 bit=0;
492 }
493 }
494 // end comment
495 if (bit) {
496 fDigitsIndexN->AddAt(i,fNDigitsN++);
497 NSidx[NNs++] =tmp;
498 }
499 }
500 }
501
502
503 if (debug) cout<<"Digits : P = "<<fNDigitsP<<" N = "<<fNDigitsN<<endl;
504
505}
506
507
508//-------------------------------------------
509
510void AliITSClusterFinderSSD::SortDigits()
511{
512
513
514//Piotr Krzysztof Skowronski
515//Warsaw University of Technology
516//skowron@if.pw.edu.pl
517
518
519
520 Int_t i;
521 if(fNDigitsP>1)
e8189707 522 for (i=0;i<fNDigitsP-1;i++)
b0f5e3fc 523 if (SortDigitsP(0,(fNDigitsP-1-i))==0) break;
524
525 if(fNDigitsN>1)
e8189707 526 for (i=0;i<fNDigitsN-1;i++)
b0f5e3fc 527 if(SortDigitsN(0,(fNDigitsN-1-i))==0) break;
528}
529
530
531
532//----------------------------------------------
533void AliITSClusterFinderSSD::FillClIndexArrays(Int_t* arrayP, Int_t *arrayN)
534{
535
536
537//Piotr Krzysztof Skowronski
538//Warsaw University of Technology
539//skowron@if.pw.edu.pl
540
541
542 register Int_t i;
e8189707 543 for (i=0; i<fNClusterP;i++)
b0f5e3fc 544 {
545 arrayP[i]=i;
546 }
e8189707 547 for (i=0; i<fNClusterN;i++)
b0f5e3fc 548 {
549 arrayN[i]=i;
550 }
551}
552
553
554//------------------------------------------------------
555void AliITSClusterFinderSSD::SortClusters(Int_t* arrayP, Int_t *arrayN)
556{
557
558
559//Piotr Krzysztof Skowronski
560//Warsaw University of Technology
561//skowron@if.pw.edu.pl
562
563
564 Int_t i;
565 if(fNClusterP>1)
e8189707 566 for (i=0;i<fNClusterP-1;i++)
b0f5e3fc 567 if (SortClustersP(0,(fNClusterP-1),arrayP)==0) break;
568
569
570 if(fNClusterN>1)
e8189707 571 for (i=0;i<fNClusterN-1;i++)
b0f5e3fc 572 if (SortClustersN(0,(fNClusterN-1),arrayN)==0) break;
573
574}
575
576
577
578//---------------------------------------------------
579Int_t AliITSClusterFinderSSD::SortClustersP(Int_t start, Int_t end, Int_t *array)
580{
581
582
583//Piotr Krzysztof Skowronski
584//Warsaw University of Technology
585//skowron@if.pw.edu.pl
586
587
588
589 Int_t right;
590 Int_t left;
591 if (start != (end - 1) ) {
592 left=this->SortClustersP(start,(start+end)/2,array);
593 right=this->SortClustersP((start+end)/2,end,array);
594 return (left || right);
595 } else {
596 left =((AliITSclusterSSD*)((*fClusterP)[array[start]]))->
597 GetDigitStripNo(0);
598 right=((AliITSclusterSSD*)((*fClusterP)[array[ end ]]))->
599 GetDigitStripNo(0);
600 if(left>right) {
601 Int_t tmp = array[start];
602 array[start]=array[end];
603 array[end]=tmp;
604 return 1;
605 } else return 0;
606 }
607
608
609}
610
611
612
613//-------------------------------------------------------
614Int_t AliITSClusterFinderSSD::SortClustersN(Int_t start, Int_t end, Int_t *array)
615{
616
617
618//Piotr Krzysztof Skowronski
619//Warsaw University of Technology
620//skowron@if.pw.edu.pl
621
622 Int_t right;
623 Int_t left;
624
625 if (start != (end - 1) ) {
626 left=this->SortClustersN(start,(start+end)/2,array);
627 right=this->SortClustersN((start+end)/2,end,array);
628 return (left || right);
629 } else {
630 left =((AliITSclusterSSD*)((*fClusterN)[array[start]]))->
631 GetDigitStripNo(0);
632 right=((AliITSclusterSSD*)((*fClusterN)[array[ end ]]))->
633 GetDigitStripNo(0);
634 if( left > right) {
635 Int_t tmp = array[start];
636 array[start]=array[end];
637 array[end]=tmp;
638 return 1;
639 } else return 0;
640 }
641
642}
643
644
645
646//-------------------------------------------------------
647void AliITSClusterFinderSSD::ClustersToPackages()
648{
649
650
651
652//Piotr Krzysztof Skowronski
653//Warsaw University of Technology
654//skowron@if.pw.edu.pl
655
656
657 Int_t *oneSclP = new Int_t[fNClusterP]; //I want to have sorted 1 S clusters
658 Int_t *oneSclN = new Int_t[fNClusterN]; //I can not sort it in TClonesArray
659 //so, I create table of indexes and
660 //sort it
661 //I do not use TArrayI on purpose
662 // MB: well, that's not true that one
663 //cannot sort objects in TClonesArray
664
e8189707 665 register Int_t i; //iterator
b0f5e3fc 666 AliITSpackageSSD *currentpkg;
667 AliITSclusterSSD *currentP;
668 AliITSclusterSSD *currentN;
669 AliITSclusterSSD *fcurrN;
670
671 Int_t lastNclIndex=0; //index of last N sidecluster
672 Int_t Nclpack=0; //number of P clusters in current package
673 Int_t Pclpack=0; //number of N clusters in current package
674 Int_t firstStripP;
675 Int_t lastStripP;
676 Int_t firstStripN;
677 Int_t lastStripN;
678 Int_t tmplastNclIndex;
679
680//Fills in One Side Clusters Index Array
681 FillClIndexArrays(oneSclP,oneSclN);
682//Sorts filled Arrays
683 SortClusters(oneSclP,oneSclN);
684
685
686 fNPackages=1;
687 new ((*fPackages)[0]) AliITSpackageSSD(fClusterP,fClusterN);
688 currentpkg = (AliITSpackageSSD*)((*fPackages)[0]);
689
690 //main loop on sorted P side clusters
e8189707 691 Int_t k;
692 for (i=0;i<fNClusterP;i++) {
b0f5e3fc 693 //Take new P side cluster
694 currentP = GetPSideCluster(oneSclP[i]);
695 currentN = GetNSideCluster(oneSclN[lastNclIndex]);
696 // take a new P cluster or not ?
697 if(IsCrossing(currentP,currentN)) {
698 // don't take a new P cluster
699 Pclpack++;
700 currentP->AddCross(oneSclN[lastNclIndex]);
701 currentN->AddCross(oneSclP[i]);
702 currentpkg->AddPSideCluster(oneSclP[i]);
703 if (Nclpack==0) {
704 currentpkg->AddNSideCluster(oneSclN[lastNclIndex]);
705 Nclpack++;
706 }
707 //check if previous N side clusters crosses with it too
708 for(k=1;k<fSFB+1;k++) {
709 //check if we are not out of array
710 if ((lastNclIndex-k)>-1) {
711 fcurrN = GetNSideCluster( oneSclN[lastNclIndex-k] );
712 if( IsCrossing(currentP,fcurrN) ) {
713 currentP->AddCross(oneSclN[lastNclIndex-k]);
714 fcurrN->AddCross(oneSclP[i]);
715 } else break; //There is no sense to check rest of clusters
716 } else break;
717 }
718 tmplastNclIndex=lastNclIndex;
719 //Check if next N side clusters crosses with it too
e8189707 720 for (Int_t k=1;k<fSFF+1;k++) {
b0f5e3fc 721 //check if we are not out of array
722 if ((tmplastNclIndex+k)<fNClusterN) {
723 fcurrN = GetNSideCluster( oneSclN[tmplastNclIndex+k] );
724 if(IsCrossing(currentP,fcurrN) ) {
725 lastNclIndex++;
726 fcurrN->AddCross(oneSclP[i]);
727 currentP->AddCross(oneSclN[tmplastNclIndex+k]);
728 currentpkg->AddNSideCluster(oneSclN[lastNclIndex]);
729 Nclpack++;
730 }else break;
731 } else break;
732 }
733
734 // end of package
735 //if( IsCrossing )
736 } else {
737 lastStripP = currentP->GetLastDigitStripNo();
738 lastStripN = currentN->GetLastDigitStripNo();
739
740 if((lastStripN-fSFF) < lastStripP ) {
741 //Take new PCluster
742 if((Nclpack>0)&& (Pclpack>0)) {
743 new ((*fPackages)[fNPackages]) AliITSpackageSSD(fClusterP,fClusterN);
744 currentpkg = (AliITSpackageSSD*)((*fPackages)[fNPackages++]);
745 }
746
747 //so we have to take another cluster on side N and check it
748 //we have to check until taken cluster's last strip will
749 //be > last of this cluster(P)
750 //so it means that there is no corresponding cluster on side N
751 //to this cluster (P)
752 //so we have to continue main loop with new "lastNclIndex"
753
754 Nclpack=0;
755 Pclpack=0;
756 //We are not sure that next N cluter will cross with this P cluster
757 //There might one, or more, clusters that do not have any
758 //corresponding cluster on the other side
759 for(;;) {
760 lastNclIndex++;
761 //Check if we are not out of array
762 if (lastNclIndex<fNClusterN) {
763 currentN = GetNSideCluster(oneSclN[lastNclIndex]);
764 if ( IsCrossing(currentP, currentN) ){
765 Nclpack++;
766 Pclpack++;
767 currentP->AddCross(oneSclN[lastNclIndex]);
768 currentN->AddCross(oneSclP[i]);
769 currentpkg->AddNSideCluster(oneSclN[lastNclIndex]);
770 currentpkg->AddPSideCluster(oneSclP[i]);
771 //Check, if next N side clusters crosses with it too
772 tmplastNclIndex=lastNclIndex;
e8189707 773 for (Int_t k=1;k<fSFF+1;k++) {
b0f5e3fc 774 //Check if we are not out of array
775 if ((tmplastNclIndex+k)<fNClusterN) {
776 fcurrN = GetNSideCluster(oneSclN[tmplastNclIndex+k]);
777 if( IsCrossing(currentP, fcurrN) ) {
778 Nclpack++;
779 lastNclIndex++;
780 currentP->AddCross(oneSclN[tmplastNclIndex+k]);
781 fcurrN->AddCross(oneSclP[i]);
782 currentpkg->AddNSideCluster(oneSclN[tmplastNclIndex+k]);
783 }else break;
784 } else break; //we are out of array
785 } // end loop
786 break;
787 } else {
788 firstStripP = currentP->GetFirstDigitStripNo();
789 firstStripN = currentN->GetFirstDigitStripNo();
790 if((firstStripN+fSFB)>=firstStripP) break;
791 else continue;
792 }
793 } else goto EndOfFunction;
794 } // end for(;;)
795 } else //EndOfPackage
796 {
797 continue;
798 } //else EndOfPackage
799
800 }//else IsCrossing
801
802 }//main for
803
804 EndOfFunction:
805 if ((Nclpack<1) || (Pclpack<1)) fNPackages--;
806
807 delete oneSclP;
808 delete oneSclN;
809
810}
811
812
813
814//-----------------------------------------------
815void AliITSClusterFinderSSD::PackagesToPoints()
816{
817
818 register Int_t i;
819 AliITSpackageSSD *currentpkg;
820 Int_t NumNcl;
821 Int_t NumPcl;
822 Int_t clusterIndex;
823 Bool_t clusterSide;
824
e8189707 825 for (i=0;i<fNPackages;i++) {
b0f5e3fc 826 //get pointer to next package
827 currentpkg = (AliITSpackageSSD*)((*fPackages)[i]);
828 NumNcl = currentpkg->GetNumOfClustersN();
829 NumPcl = currentpkg->GetNumOfClustersP();
830 if(debug) cout<<"New Package\nNumPcl ="<<NumPcl<<" NumNcl ="<<NumNcl<<"\n";
831
832 while(((NumPcl>=2)&&(NumNcl>2))||((NumPcl>2)&&(NumNcl>=2))) {
833 //This is case of "big" pacakage
834 //conditions see 3 lines above
835 //if in the big package exists cluster with one cross
836 //we can reconstruct this point without any geometrical ambiguities
837 if(currentpkg->GetClusterWithOneCross(clusterIndex, clusterSide) )
838 {
839 ResolveClusterWithOneCross(currentpkg, clusterIndex, clusterSide);
840 } else if (clusterIndex == -2) {
841 NumPcl = 0;
842 NumNcl = 0;
843 break;
844 } else {
845 if ( (NumNcl==NumPcl) && (NumPcl<10)) {
846 //if there is no cluster with one cross
847 //we can resolve whole package at once
848 //by finding best combination
849 //we can make combinations when NumNcl==NumPcl
850 if (ResolvePackageBestCombin(currentpkg)) {
851 //sometimes creating combinations fail,
852 //but it happens very rarely
853 NumPcl = 0;
854 NumNcl = 0;
855 break;
856 } else ResolveOneBestMatchingPoint(currentpkg);
857 } else {
858 ResolveOneBestMatchingPoint(currentpkg);
859 }
860 }
861 NumNcl = currentpkg->GetNumOfClustersN();
862 NumPcl = currentpkg->GetNumOfClustersP();
863 } // end while
864 if ((NumPcl==1)&&(NumNcl==1)) {
865 ResolveSimplePackage(currentpkg);
866 continue;
867 }
868 if (NumPcl==1) {
869 ResolvePackageWithOnePSideCluster(currentpkg);
870 continue;
871 }
872 if (NumNcl==1) {
873 ResolvePackageWithOneNSideCluster(currentpkg);
874 continue;
875 }
876 if ((NumPcl==2)&&(NumNcl==2)) {
877 ResolveTwoForTwoPackage(currentpkg);
878 continue;
879 }
880 if ((NumPcl==0)&&(NumNcl>0)) { }
881 if ((NumNcl==0)&&(NumPcl>0)) { }
882
883 } // end loop over fNPackages
884
885
886}
887
888
889//----------------------------------------------------------
890
891void AliITSClusterFinderSSD::
892ResolveClusterWithOneCross(AliITSpackageSSD *currentpkg, Int_t clusterIndex, Bool_t clSide)
893{
894
e8189707 895 if (clSide == fgkSIDEP) ResolvePClusterWithOneCross(currentpkg,clusterIndex);
b0f5e3fc 896 else ResolveNClusterWithOneCross(currentpkg,clusterIndex);
897
898}
899
900
901//---------------------------------------------------------
902void AliITSClusterFinderSSD::
903ResolvePClusterWithOneCross(AliITSpackageSSD *pkg, Int_t clusterIndex)
904{
905
906/*
907There is cluster (side P) which crosses with only one cluster on the other
908side (N)
909
910ie:
911 \ / \/ /
912 \ / /\ /
913 \/ / \/
914 /\ / /\
915 / \/ / \ .....
916 / /\ / \
917 / / \/ \
918 / / /\ \
919 / / / \ \
920*/
921
922
923 AliITSclusterSSD * clusterP;
924 AliITSclusterSSD * clusterN;
925
926 Float_t posClusterP; //Cluster P position in strip coordinates
927 Float_t posClusterN; //Cluster N position in strip coordinates
928
929 Float_t posErrorClusterP;
930 Float_t posErrorClusterN;
931
932 Float_t sigClusterP;
933 Float_t sigClusterN;
934
935 Float_t sigErrorClusterP;
936 Float_t sigErrorClusterN;
937
938 Float_t ps;
939 Float_t ns;
940
941 Float_t Chicomb;
942 Int_t clusterIdx;
943
944 if (debug) cout<<"ResolvePClusterWithOneCross\n";
945
946 clusterP=pkg->GetPSideCluster(clusterIndex);
947 posClusterP = GetClusterZ(clusterP);
948 posErrorClusterP = clusterP->GetPositionError();
949 sigClusterP = ps = clusterP->GetTotalSignal();
950 sigErrorClusterP = clusterP->GetTotalSignalError();
951 //carefully, it returns index in TClonesArray
952
953 clusterIdx = clusterP->GetCross(0);
954 clusterN=this->GetNSideCluster(clusterIdx);
955 ns = clusterN->GetTotalSignal();
956 posClusterN = GetClusterZ(clusterN);
957 posErrorClusterN = clusterN->GetPositionError();
e8189707 958 pkg->DelCluster(clusterIndex,fgkSIDEP);
be33dccb 959 sigClusterN = ps/fPNsignalRatio;
b0f5e3fc 960 // there is no sonse to check how signal ratio is far from perfect
961 // matching line if the if below it is true
962 if (ns < sigClusterN) {
963 sigClusterN=ns;
be33dccb 964 if (debug) cout<<"n1 < p1/fPNsignalRatio";
b0f5e3fc 965 if (debug) cout<<"Attempting to del cluster N "<<clusterIdx<<" ... \n";
e8189707 966 pkg->DelClusterOI(clusterIdx,fgkSIDEN);
b0f5e3fc 967 } else {
968 //Let's see how signal ratio is far from perfect matching line
969 Chicomb = DistToPML(ps,ns);
970 if (debug) cout<<"Chic "<<Chicomb<<"\n";
971 if (Chicomb > falpha2) {
972 //it is near, so we can risk throwing this cluster away too
973 if (debug) cout<<"Attempting to del cluster N "<<clusterIdx<<"...\n";
e8189707 974 pkg->DelClusterOI(clusterIdx,fgkSIDEN);
b0f5e3fc 975 } else {
976 clusterN->CutTotalSignal(sigClusterN);
977 if (debug) cout <<"Signal cut |||||||||||||\n";
978 }
979 }
980 sigErrorClusterN= clusterN->GetTotalSignalError();
981 CreateNewRecPoint(posClusterP,posErrorClusterP,posClusterN,posErrorClusterN,
982 sigClusterP+sigClusterN, sigErrorClusterN+sigErrorClusterP,
983 clusterP, clusterN, 0.75);
984
985}
986
987
988//------------------------------------------------------
989void AliITSClusterFinderSSD::
990ResolveNClusterWithOneCross(AliITSpackageSSD *pkg, Int_t clusterIndex)
991{
992
993 AliITSclusterSSD * clusterP;
994 AliITSclusterSSD * clusterN;
995
996 Float_t posClusterP; //Cluster P position in strip coordinates
997 Float_t posClusterN; //Cluster N position in strip coordinates
998
999 Float_t posErrorClusterP;
1000 Float_t posErrorClusterN;
1001
1002 Float_t sigClusterP;
1003 Float_t sigClusterN;
1004
1005 Float_t sigErrorClusterP;
1006 Float_t sigErrorClusterN;
1007
1008 Float_t ps;
1009 Float_t ns;
1010
1011 Float_t Chicomb;
1012 Int_t clusterIdx;
1013
1014 if (debug) cout<<"ResolveNClusterWithOneCross\n";
1015
1016 clusterN=pkg->GetNSideCluster(clusterIndex);
1017 posClusterN = GetClusterZ(clusterN);
1018 posErrorClusterN = clusterN->GetPositionError();
1019 sigClusterN = ns = clusterN->GetTotalSignal();
1020 sigErrorClusterN = clusterN->GetTotalSignalError();
1021 //carefully, it returns index in TClonesArray
1022 clusterIdx = clusterN->GetCross(0);
1023 clusterP=this->GetPSideCluster(clusterIdx);
1024 ps = clusterP->GetTotalSignal();
1025 posClusterP = GetClusterZ(clusterP);
1026 posErrorClusterP = clusterP->GetPositionError();
e8189707 1027 pkg->DelCluster(clusterIndex,fgkSIDEN);
be33dccb 1028 sigClusterP=ns*fPNsignalRatio;
b0f5e3fc 1029 // there is no sonse to check how signal ratio is far from perfect
1030 // matching line if the if below it is true
1031 if (ps < sigClusterP) {
1032 sigClusterP = ps;
be33dccb 1033 if (debug) cout<<"ps < ns*fPNsignalRatio";
b0f5e3fc 1034 if (debug) cout<<"Attempting to del cluster P "<<clusterIdx<<" ... \n";
e8189707 1035 pkg->DelClusterOI(clusterIdx,fgkSIDEP);
b0f5e3fc 1036 } else {
1037 //Let's see how signal ratio is far from perfect matching line
1038 Chicomb = DistToPML(ps,ns);
1039 if (debug) cout<<"Chic "<<Chicomb<<"\n";
1040 if (Chicomb > falpha2) {
1041 //it is near, so we can risk frowing this cluster away too
1042 if (debug) cout<<"Attempting to del cluster P "<<clusterIdx<<"...\n";
e8189707 1043 pkg->DelClusterOI(clusterIdx,fgkSIDEP);
b0f5e3fc 1044 } else {
1045 clusterN->CutTotalSignal(sigClusterP);
1046 if (debug) cout <<"Signal cut ------------\n";
1047 }
1048 }
1049 sigErrorClusterP= clusterP->GetTotalSignalError();
1050 CreateNewRecPoint( posClusterP,posErrorClusterP,posClusterN,posErrorClusterN,
1051 sigClusterP+sigClusterN,sigErrorClusterN+sigErrorClusterP,
1052 clusterP, clusterN, 0.75);
1053
1054
1055}
1056
1057
1058
1059//-------------------------------------------------------
1060
1061Bool_t AliITSClusterFinderSSD::
1062ResolvePackageBestCombin(AliITSpackageSSD *pkg)
1063{
1064 if (debug) cout<<"NumNcl==NumPcl ("<<pkg->GetNumOfClustersN()
1065 <<"=="<<pkg->GetNumOfClustersP()<<"); Generating combinations ... \n";
1066
1067
1068 AliITSclusterSSD * clusterP;
1069 AliITSclusterSSD * clusterN;
1070
1071 Int_t Ncombin; //Number of combinations
1072 Int_t itera; //iterator
1073 Int_t sizet=1; //size of array to allocate
e8189707 1074
b0f5e3fc 1075 Int_t NP = pkg->GetNumOfClustersP();
e8189707 1076 for (itera =2; itera <= NP ;itera ++) {
b0f5e3fc 1077 sizet=sizet*itera;
1078 if (sizet > 10000) {
1079 sizet=10000;
1080 break;
1081 }
1082 }
1083
e8189707 1084 Int_t** combin = new Int_t*[sizet]; //2D array to keep combinations in
b0f5e3fc 1085
e8189707 1086 for (itera =0; itera <sizet;itera++) {
b0f5e3fc 1087 combin[itera] = new Int_t[NP+1];
1088 }
1089
1090 pkg->GetAllCombinations(combin,Ncombin,sizet);
1091
1092 if (Ncombin==0) {
1093 if (debug) cout<<"No combin Error";
1094 return kFALSE;
1095 }
1096//calculate which combination fits the best to perfect matching line
1097 Int_t bc = GetBestComb(combin,Ncombin,NP,pkg);
1098 if (debug) cout<<" bc "<<bc <<"\n";
1099
e8189707 1100 for (itera =0; itera < NP; itera ++) {
b0f5e3fc 1101 clusterP = pkg->GetPSideCluster(itera);
1102 //carefully here
1103 //becase AliITSclusterSSD::GetCross returns index in
1104 //AliITSmoduleSSD.fClusterP, which is put to "combin"
1105 clusterN = GetNSideCluster(combin[bc][itera]);
1106 CreateNewRecPoint(clusterP, clusterN, 0.75);
1107 }
1108
e8189707 1109 for (itera =0; itera <sizet;itera++) {
b0f5e3fc 1110 delete [](combin[itera]);
1111 }
1112 delete [] combin;
1113 return kTRUE;
1114
1115}
1116
1117
1118//----------------------------------------------------
1119void AliITSClusterFinderSSD::
1120ResolveOneBestMatchingPoint(AliITSpackageSSD *pkg)
1121{
1122
1123 Int_t ni, pi;
1124
1125 Int_t prvP, prvN;
1126 Int_t nextP, nextN=0;
1127
1128
1129 AliITSclusterSSD * clusterP;
1130 AliITSclusterSSD * clusterN;
1131
1132 Bool_t split = kFALSE;
1133
1134 if (debug) cout<<"ResolveOneBestMatchingPoint \n";
1135
1136 GetBestMatchingPoint(pi, ni, pkg);
1137
1138 clusterP = GetPSideCluster(pi);
1139 clusterN = GetNSideCluster(ni);
1140
1141 CreateNewRecPoint(clusterP, clusterN, 0.75);
1142
1143 if ((nextP=pkg->GetNextPIdx(pi))!=-1)
1144 if ((nextN=pkg->GetNextNIdx(ni))!=-1)
1145 if ((prvP=pkg->GetPrvPIdx(pi))!=-1)
1146 if ((prvN=pkg->GetPrvNIdx(ni))!=-1)
1147 if( !(GetPSideCluster(prvP)->IsCrossingWith(nextN)))
1148 if( !(GetPSideCluster(nextP)->IsCrossingWith(prvN)))
1149 {
1150 split = kTRUE;
1151 }
1152
1153
e8189707 1154 pkg->DelClusterOI(pi, fgkSIDEP);
1155 pkg->DelClusterOI(ni, fgkSIDEN);
b0f5e3fc 1156
1157 if (split) {
1158 if (debug) cout<<"spltting package ...\n";
1159 new ((*fPackages)[fNPackages]) AliITSpackageSSD(fClusterP,fClusterN);
1160 pkg->
1161 SplitPackage(nextP,nextN,(AliITSpackageSSD*)((*fPackages)[fNPackages++]));
1162 }
1163
1164}
1165
1166
1167//--------------------------------------------------
1168void AliITSClusterFinderSSD::ResolveSimplePackage(AliITSpackageSSD *pkg)
1169{
1170
1171 AliITSclusterSSD * clusterP;
1172 AliITSclusterSSD * clusterN;
1173
1174 clusterP = pkg->GetPSideCluster(0);
1175
1176 clusterN = pkg->GetNSideCluster(0);
1177
1178 CreateNewRecPoint(clusterP, clusterN, 1.0);
1179
1180
1181}
1182
1183
1184//--------------------------------------------------
1185void AliITSClusterFinderSSD:: ResolvePackageWithOnePSideCluster(AliITSpackageSSD *pkg) {
1186
1187
1188/*
1189 \ \ \ /
1190 \ \ \/
1191 \ \ /\
1192 \ \/ \
1193 \ /\ \
1194 \/ \ \
1195 /\ \ \
1196 / \ \ \
1197
1198*/
1199
1200
1201/************************/
1202/**************************/
1203/*XP SP itg jest tylko jeden nie musi byc tablica i w peetli nie trzeba po niej iterowcac*/
1204/***************************/
1205
1206
1207 Int_t k;
1208 AliITSclusterSSD * clusterP;
1209 AliITSclusterSSD * clusterN;
1210
1211 Int_t NN=pkg->GetNumOfClustersN();
1212 Float_t sumsig = 0;
1213
1214 Float_t XP;
1215 Float_t XPerr;
1216
1217 Float_t * XN = new Float_t[NN];
1218 Float_t * XNerr = new Float_t[NN];
1219
1220 Float_t * SP = new Float_t[NN];
1221 Float_t * SN = new Float_t[NN];
1222
1223 Float_t * SPerr = new Float_t[NN];
1224 Float_t * SNerr = new Float_t[NN];
1225
1226 Float_t p1;
1227 Float_t p1err;
1228
1229 clusterP = pkg->GetPSideCluster(0);
1230 p1 = clusterP->GetTotalSignal();
1231
1232 XP = GetClusterZ(clusterP);
1233 XPerr = clusterP->GetPositionError();
1234 p1err = clusterP->GetTotalSignalError();
1235
e8189707 1236 for (k=0;k<NN;k++) {
b0f5e3fc 1237 SN[k] = pkg->GetNSideCluster(k)->GetTotalSignal();
1238 SNerr[k] = pkg->GetNSideCluster(k)->GetTotalSignalError();
1239 sumsig += SN[k];
1240 }
e8189707 1241 for (k=0;k<NN;k++) {
b0f5e3fc 1242 clusterN = pkg->GetNSideCluster(k);
1243 SP[k]= p1*SN[k]/sumsig;
1244 SPerr[k] = p1err*SN[k]/sumsig;
1245 XN[k]=GetClusterZ(clusterN);
1246 XNerr[k]= clusterN->GetPositionError();
1247 CreateNewRecPoint(XP,XPerr, XN[k],XNerr[k], SP[k]+SN[k], SPerr[k]+SNerr[k], clusterP, clusterN, 1.0);
1248
1249 }
1250
1251}
1252
1253
1254//---------------------------------------------------------
1255void AliITSClusterFinderSSD::ResolvePackageWithOneNSideCluster(AliITSpackageSSD *pkg) {
1256
1257
1258/*
1259 \ / / /
1260 \ / / /
1261 \/ / /
1262 /\ / /
1263 / \/ /
1264 / /\ /
1265 / / \/
1266 / / /\
1267 / / / \
1268*/
1269
1270
1271 Int_t k;
1272 AliITSclusterSSD * clusterP;
1273 AliITSclusterSSD * clusterN;
1274
1275 Int_t NP=pkg->GetNumOfClustersP();
1276 Float_t sumsig = 0;
1277
1278 Float_t * XP = new Float_t[NP];
1279 Float_t * XPerr = new Float_t[NP];
1280
1281 Float_t XN;
1282 Float_t XNerr;
1283
1284 Float_t * SP = new Float_t[NP];
1285 Float_t * SN = new Float_t[NP];
1286
1287 Float_t * SPerr = new Float_t[NP];
1288 Float_t * SNerr = new Float_t[NP];
1289
1290 Float_t n1;
1291 Float_t n1err;
1292
1293 clusterN = pkg->GetNSideCluster(0);
1294 n1 = clusterN->GetTotalSignal();
1295
1296 XN = GetClusterZ(clusterN);
1297 XNerr = clusterN->GetPositionError();
1298
1299 n1err=clusterN->GetTotalSignalError();
1300
e8189707 1301 for (k=0;k<NP;k++) {
b0f5e3fc 1302 SP[k] = pkg->GetPSideCluster(k)->GetTotalSignal();
1303 sumsig += SP[k];
1304 SPerr[k] = pkg->GetPSideCluster(k)->GetTotalSignalError();
1305 }
1306
e8189707 1307 for (k=0;k<NP;k++) {
b0f5e3fc 1308 clusterP = pkg->GetPSideCluster(k);
1309 SN[k]= n1*SP[k]/sumsig;
1310 XP[k]=GetClusterZ(clusterP);
1311 XPerr[k]= clusterP->GetPositionError();
1312 SNerr[k] = n1err*SP[k]/sumsig;
1313 CreateNewRecPoint(XP[k],XPerr[k], XN,XNerr, SP[k]+SN[k], SPerr[k]+SNerr[k],clusterP, clusterN, 1.0);
1314 }
1315
1316
1317}
1318
1319/**********************************************************************/
1320/********* 2 X 2 *********************************************/
1321/**********************************************************************/
1322
1323void AliITSClusterFinderSSD::
1324ResolveTwoForTwoPackage(AliITSpackageSSD *pkg)
1325{
1326
1327 AliITSclusterSSD *clusterP1 = pkg->GetPSideCluster(0);
1328 AliITSclusterSSD *clusterP2 = pkg->GetPSideCluster(1);
1329 AliITSclusterSSD *clusterN1 = pkg->GetNSideCluster(0);
1330 AliITSclusterSSD *clusterN2 = pkg->GetNSideCluster(1);
1331
1332 AliITSclusterSSD *tmp;
1333
1334 Float_t p1sig;
1335 Float_t p2sig;
1336 Float_t n1sig;
1337 Float_t n2sig;
1338
1339 Float_t p1sigErr;
1340 Float_t p2sigErr;
1341 Float_t n1sigErr;
1342 Float_t n2sigErr;
1343
1344 Float_t ZposP1;
1345 Float_t ZposP2;
1346 Float_t ZposN1;
1347 Float_t ZposN2;
1348
1349 Float_t ZposErrP1;
1350 Float_t ZposErrP2;
1351 Float_t ZposErrN1;
1352 Float_t ZposErrN2;
1353
1354 Float_t D12;
1355
1356 Float_t Chicomb1;
1357 Float_t Chicomb2;
1358
1359 if(clusterP1->GetDigitStripNo(0) > clusterP2->GetDigitStripNo(0)) {
1360 if (debug) cout<<"P strips flip\n";
1361 tmp = clusterP1;
1362 clusterP1 = clusterP2;
1363 clusterP2 = tmp;
1364 }
1365 if(clusterN1->GetDigitStripNo(0) > clusterN2->GetDigitStripNo(0)) {
1366 if (debug) cout<<"N strips flip\n";
1367 tmp = clusterN1;
1368 clusterN1 = clusterN2;
1369 clusterN2 = tmp;
1370 }
1371 p1sig = clusterP1->GetTotalSignal();
1372 p2sig = clusterP2->GetTotalSignal();
1373 n1sig = clusterN1->GetTotalSignal();
1374 n2sig = clusterN2->GetTotalSignal();
1375
1376
1377 p1sigErr = clusterP1->GetTotalSignalError();
1378 n1sigErr = clusterN1->GetTotalSignalError();
1379 p2sigErr = clusterP2->GetTotalSignalError();
1380 n2sigErr = clusterN2->GetTotalSignalError();
1381
1382 ZposP1 = clusterP1->GetPosition();
1383 ZposN1 = clusterN1->GetPosition();
1384 ZposP2 = clusterP2->GetPosition();
1385 ZposN2 = clusterN2->GetPosition();
1386
1387 ZposErrP1 = clusterP1->GetPositionError();
1388 ZposErrN1 = clusterN1->GetPositionError();
1389 ZposErrP2 = clusterP2->GetPositionError();
1390 ZposErrN2 = clusterN2->GetPositionError();
1391
1392 //in this may be two types:
1393 // 1.When each cluster crosses with 2 clusters on the other side
1394 // gives 2 ghosts and 2 real points
1395 //
1396 // 2.When two clusters (one an each side) crosses with only one on
1397 // the other side and two crosses (one on the each side) with
1398 // two on the other gives 2 or 3 points,
1399
1400 if (debug) cout<<"2 for 2 ambiguity ...";
1401 /***************************/
1402 /*First sort of combination*/
1403 /***************************/
1404
1405 if((clusterP1->GetCrossNo()==2) && (clusterN1->GetCrossNo()==2)) {
1406 if (debug) cout<<"All clusters has 2 crosses\n";
1407 D12 = TMath::Sqrt((Float_t)((p1sig -p2sig)*(p1sig -p2sig) + (n1sig -n2sig)*(n1sig -n2sig)));
1408
1409 if(debug) cout<<"Distance between points in P(N) plane D12 = "<<D12<<"\n";
1410
1411 /*********************************************/
1412 /*resolving ambiguities: */
1413 /*we count Chicomb's and we take combination */
1414 /*giving lower Chicomb as a real points */
1415 /*Keep only better combinantion */
1416 /*********************************************/
1417
1418 if (D12 > (falpha3*17.8768)) {
1419 if (debug) cout<<"decided to take only one pair \n";
1420 Chicomb1 = DistToPML(p1sig,n1sig) + DistToPML(p2sig,n2sig);
1421 Chicomb2 = DistToPML(p2sig,n1sig) + DistToPML(p1sig,n2sig);
1422 if (debug) {
1423 cout<<" 00 11 combination : "<<Chicomb1<<" 01 10 combination : "<<Chicomb2<<" \n";
1424 cout<<"p1 = "<<p1sig<<" n1 = "<<n1sig<<" p2 = "<<p2sig<<" n2 = "<<n2sig<<"\n";
1425 }
1426 if (Chicomb1 < Chicomb2) {
1427 if (debug) cout<<"00 11";
1428 CreateNewRecPoint(ZposP1,ZposErrP1, ZposN1,ZposErrN1, p1sig+n1sig, p1sigErr+n1sigErr, clusterP1, clusterN1, 0.75);
1429 CreateNewRecPoint(ZposP2,ZposErrP2, ZposN2,ZposErrN2, p2sig+n2sig, p2sigErr+n2sigErr, clusterP2, clusterN2, 0.75);
1430
1431
1432 //second cominantion
1433 } else {
1434 if (debug) cout<<"01 10";
1435 CreateNewRecPoint(ZposP1,0, ZposN2,0, p1sig+n2sig, p1sigErr+n2sigErr, clusterP1, clusterN2, 0.75);
1436 CreateNewRecPoint(ZposP2,0, ZposN1,0, p2sig+n1sig, p2sigErr+n1sigErr, clusterP2, clusterN1, 0.75);
1437 } //end second combinantion
1438 //if (D12 > falpha3*17.8768)
1439 //keep all combinations
1440 } else {
1441 if (debug) cout<<"We decide to take all points\n";
1442 CreateNewRecPoint(ZposP1,ZposErrP1, ZposN1,ZposErrN1, p1sig+n1sig, p1sigErr+n1sigErr, clusterP1, clusterN1, 0.5);
1443 CreateNewRecPoint(ZposP2,ZposErrP2, ZposN2,ZposErrN2, p2sig+n2sig, p2sigErr+n2sigErr, clusterP2, clusterN2, 0.5);
1444 CreateNewRecPoint(ZposP1,ZposErrP1, ZposN2,ZposErrN2, p1sig+n2sig, p1sigErr+n2sigErr, clusterP1, clusterN2, 0.5);
1445 CreateNewRecPoint(ZposP2,ZposErrP2, ZposN1,ZposErrN1, p2sig+n1sig, p2sigErr+n1sigErr, clusterP2, clusterN1, 0.5);
1446 }
1447 }
1448 // ad.2 Second type of combination
1449 else {
1450 Chicomb1 = DistToPML(p1sig,n1sig) + DistToPML(p2sig,n2sig);
1451 if (debug) cout<<"\nhere can be reconstructed 3 points: chicomb = "<<Chicomb1<<"\n";
1452
1453 if (Chicomb1<falpha1) {
1454 if (debug) cout<<"\nWe decided to take 3rd point";
1455 if (clusterP1->GetCrossNo()==1) {
1456 if (debug) cout<<"... P1 has one cross\n";
be33dccb 1457 n1sig = p1sig/fPNsignalRatio;
1458 p2sig = n2sig*fPNsignalRatio;
b0f5e3fc 1459
1460 clusterN1->CutTotalSignal(n1sig);
1461 clusterP2->CutTotalSignal(p2sig);
1462
1463 CreateNewRecPoint(ZposP2,ZposErrP2, ZposN1,ZposErrN1, p2sig+n1sig, p2sigErr+n1sigErr, clusterP2, clusterN1, 0.5);
1464
1465 n1sig = clusterN1->GetTotalSignal();
1466 p2sig = clusterP2->GetTotalSignal();
1467
1468 } else {
1469 if (debug) cout<<"... N1 has one cross\n";
1470
be33dccb 1471 n2sig=p2sig/fPNsignalRatio;
1472 p1sig=n1sig*fPNsignalRatio;
b0f5e3fc 1473
1474 clusterN2->CutTotalSignal(n2sig);
1475 clusterP1->CutTotalSignal(p1sig);
1476
1477 CreateNewRecPoint(ZposP1,ZposErrP1, ZposN2,ZposErrN2, p1sig+n2sig, p1sigErr+n2sigErr, clusterP1, clusterN2, 0.5);
1478
1479 n2sig=clusterN2->GetTotalSignal();
1480 p1sig=clusterP1->GetTotalSignal();
1481 }
1482 } else {
1483 if (debug) cout<<"\nWe decided NOT to take 3rd point\n";
1484 }
1485
1486 CreateNewRecPoint(ZposP1,ZposErrP1, ZposN1,ZposErrN1, p1sig+n1sig, p1sigErr+n1sigErr, clusterP1, clusterN1, 1.0);
1487 CreateNewRecPoint(ZposP2,ZposErrP2, ZposN2,ZposErrN2, p2sig+n2sig, p2sigErr+n2sigErr, clusterP2, clusterN2, 1.0);
1488
1489 }
1490
1491}
1492
1493
1494
1495//------------------------------------------------------
1496Bool_t AliITSClusterFinderSSD::
1497CreateNewRecPoint(Float_t P, Float_t dP, Float_t N, Float_t dN,
1498 Float_t Sig,Float_t dSig,
1499 AliITSclusterSSD *clusterP, AliITSclusterSSD *clusterN,
1500 Stat_t prob)
1501{
1502
1503 const Float_t kdEdXtoQ = 2.778e+8;
1504 const Float_t kconv = 1.0e-4;
1505 const Float_t kRMSx = 20.0*kconv; // microns->cm ITS TDR Table 1.3
1506 const Float_t kRMSz = 830.0*kconv; // microns->cm ITS TDR Table 1.3
1507
1508
1509 Float_t p=P;
1510 Float_t n=N;
e8189707 1511 Int_t stripP, stripN;
1512 Int_t sigP, sigN;
1513 AliITSdigitSSD *digP, *digN;
b0f5e3fc 1514 if (GetCrossing(P,N)) {
1515 GetCrossingError(dP,dN);
1516 AliITSRawClusterSSD cnew;
1517 Int_t nstripsP=clusterP->GetNumOfDigits();
1518 Int_t nstripsN=clusterN->GetNumOfDigits();
1519 cnew.fMultiplicity=nstripsP;
1520 cnew.fMultiplicityN=nstripsN;
1521 /*
1522 if (nstripsP > 100) {
1523 printf("multiplicity > 100 - increase the dimension of the arrays %d\n",nstripsP);
1524 nstripsP=100;
1525 }
e8189707 1526 for(int i=0;i<nstripsP;i++) {
b0f5e3fc 1527 // check if 'clusterP->GetDigitStripNo(i)' returns the digit index
1528 cnew.fIndexMap[i] = clusterP->GetDigitStripNo(i);
1529 }
1530 if (nstripsN > 100) {
1531 printf("multiplicity > 100 - increase the dimension of the arrays %d\n",nstripsN);
1532 nstripsN=100;
1533 }
e8189707 1534 for(int i=0;i<nstripsN;i++) {
b0f5e3fc 1535 // check if 'clusterN->GetDigitStripNo(i)' returns the digit index
1536 cnew.fIndexMapN[i] = clusterN->GetDigitStripNo(i);
1537 }
1538 */
1539 cnew.fQErr=dSig;
1540 //cnew.fProbability=(float)prob;
1541 fITS->AddCluster(2,&cnew);
e8189707 1542 fSegmentation->GetPadIxz(P,N,stripP,stripN);
1543 digP = (AliITSdigitSSD*)fMap->GetHit(1,stripP);
1544 digN = (AliITSdigitSSD*)fMap->GetHit(0,stripN);
1545 sigP = digP->fSignal;
1546 sigN = digN->fSignal;
b0f5e3fc 1547 // add the rec point info
1548 AliITSRecPoint rnew;
1549 rnew.SetX(P*kconv);
1550 rnew.SetZ(N*kconv);
1551 rnew.SetQ(Sig);
1552 rnew.SetdEdX(Sig/kdEdXtoQ);
1553 rnew.SetSigmaX2( kRMSx* kRMSx);
1554 rnew.SetSigmaZ2( kRMSz* kRMSz);
1555 rnew.SetProbability((float)prob);
e8189707 1556 if(sigP > sigN) {
1557 rnew.fTracks[0]=digP->fTracks[0];
1558 rnew.fTracks[1]=digP->fTracks[1];
1559 rnew.fTracks[2]=digP->fTracks[2];
1560 } else {
1561 rnew.fTracks[0]=digN->fTracks[0];
1562 rnew.fTracks[1]=digN->fTracks[1];
1563 rnew.fTracks[2]=digN->fTracks[2];
1564 }
b0f5e3fc 1565 fITS->AddRecPoint(rnew);
1566 /*
1567 // it was
1568 fPointsM->AddLast( (TObject*)
1569 new AliITSRecPointSSD(P,dP,0,0,N,dN,Sig,dSig,fLayer,fLad,fDet,clusterP,clusterN,prob) );
1570 */
1571
1572 if(debug) cout<<"\n"<<": ImpactPoint Created: X = "<<P<<" Z = "<<N<<"; P = "<<p<<" N = "<<n<<"\n";
1573 return kTRUE;
1574 } else {
1575 if (debug) {
1576 cout<<"\n"<<": ATENTION : stupid ImpactPoit Point X = "<<P<<" Z = "<<N<<"; P = "<<p<<" N = "<<n<<"\n";
1577 }
1578 return kFALSE;
1579 }
1580}
1581
1582
1583
1584/**********************************************************************/
1585Bool_t AliITSClusterFinderSSD::
1586CreateNewRecPoint(AliITSclusterSSD *clusterP, AliITSclusterSSD *clusterN, Stat_t prob)
1587{
1588 Float_t posClusterP; //Cluster P position in strip coordinates
1589 Float_t posClusterN; //Cluster N position in strip coordinates
1590
1591 Float_t posErrorClusterP;
1592 Float_t posErrorClusterN;
1593
1594 Float_t sigClusterP;
1595 Float_t sigClusterN;
1596
1597 Float_t sigErrorClusterP;
1598 Float_t sigErrorClusterN;
1599
1600 posClusterP = clusterP->GetPosition();
1601 posErrorClusterP = clusterP->GetPositionError();
1602 sigClusterP = clusterP->GetTotalSignal();
1603 sigErrorClusterP = clusterP->GetTotalSignalError();
1604
1605 posClusterN = clusterN->GetPosition();
1606 posErrorClusterN = clusterN->GetPositionError();
1607 sigClusterN = clusterN->GetTotalSignal();
1608 sigErrorClusterN = clusterN->GetTotalSignalError();
1609
1610 return CreateNewRecPoint( posClusterP,posErrorClusterP, posClusterN,posErrorClusterN,
1611 sigClusterP+sigClusterN, sigErrorClusterN+sigErrorClusterP,
1612 clusterP, clusterN, prob);
1613
1614}
1615
1616//--------------------------------------------------
1617Bool_t AliITSClusterFinderSSD::IsCrossing(AliITSclusterSSD* p, AliITSclusterSSD* n)
1618{
1619 Float_t x = p->GetPosition();
1620 Float_t y = n->GetPosition();
1621 return GetCrossing(x,y);
1622}
1623
1624
1625//----------------------------------------------
1626Bool_t AliITSClusterFinderSSD::Strip2Local( Float_t stripP, Float_t stripN, Float_t &Z,Float_t &X)
1627{
1628
1629
1630//Piotr Krzysztof Skowronski
1631//Warsaw University of Technology
1632//skowron@if.pw.edu.pl
1633
1634/*
1635 Z = (stripN-stripP)*fFactorOne;
1636 X = (stripN + stripP - fStripZeroOffset)*fFactorTwo;
1637*/
1638
1639 Float_t P =stripP;
1640 Float_t N =stripN;
1641
1642 GetCrossing(P,N);
1643 X=P;
1644 Z=N;
1645 if (debug) cout<<"P="<<stripP<<" N="<<stripN<<" X = "<<X<<" Z = "<<Z<<"\n";
e8189707 1646 if ((Z<2.1)&&(Z>-2.1)&&(X<3.65)&&(X>-3.65)) return kTRUE;
1647 else return kFALSE;
b0f5e3fc 1648
1649}
1650
1651
1652//-----------------------------------------------------------
1653Float_t AliITSClusterFinderSSD::GetClusterZ(AliITSclusterSSD* clust)
1654{
1655
1656 return clust->GetPosition();
1657
1658}
1659
1660//-------------------------------------------------------------
1661Int_t AliITSClusterFinderSSD::GetBestComb
1662(Int_t** comb,Int_t Ncomb, Int_t Ncl, AliITSpackageSSD * pkg)
1663{
1664
1665
1666//Piotr Krzysztof Skowronski
1667//Warsaw University of Technology
1668//skowron@if.pw.edu.pl
1669
1670//returns index of best combination in "comb"
1671//comb : sets of combinations
1672// in given combination on place "n" is index of
1673// Nside cluster coresponding to "n" P side cluster
1674//
1675//Ncomb : number of combinations == number of rows in "comb"
1676//
1677//Ncl : number of clusters in each row == number of columns in "comb"
1678//
1679//pkg : package
1680
1681
1682 Float_t currbval=-1; //current best value,
1683 //starting value is set to number negative,
1684 //to be changed in first loop turn automatically
1685
1686 Int_t curridx=0; //index of combination giving currently the best chi^2
1687 Float_t chi;
1688 Float_t ps, ns; //signal of P cluster and N cluster
1689
e8189707 1690 for (Int_t i=0;i<Ncomb;i++)
b0f5e3fc 1691 {
1692 chi=0;
1693
e8189707 1694 for (Int_t j=0;j<Ncl;j++)
b0f5e3fc 1695 {
1696 ps = pkg->GetPSideCluster(j)->GetTotalSignal(); //carrefully here, different functions
1697 ns = GetNSideCluster(comb[i][j])->GetTotalSignal();
1698 chi+=DistToPML(ps,ns);
1699 }
1700
1701 if (currbval==-1) currbval = chi;
1702 if (chi<currbval)
1703 {
1704 curridx = i;
1705 currbval =chi;
1706 }
1707 }
1708
1709
1710 return curridx;
1711
1712}
1713
1714//--------------------------------------------------
1715void AliITSClusterFinderSSD::GetBestMatchingPoint
1716(Int_t & ip, Int_t & in, AliITSpackageSSD* pkg )
1717{
1718
1719
1720//Piotr Krzysztof Skowronski
1721//Warsaw University of Technology
1722//skowron@if.pw.edu.pl
1723
1724 if (debug) pkg->PrintClusters();
1725
1726 Float_t ps, ns; //signals on side p & n
1727
1728 Int_t nc; // number of crosses in given p side cluster
1729 Int_t n,p;
1730 AliITSclusterSSD *curPcl, *curNcl; //current p side cluster
1731 Float_t bestchi, chi;
1732 ip=p=pkg->GetPSideClusterIdx(0);
1733 in=n=pkg->GetNSideClusterIdx(0);
1734
1735 bestchi=DistToPML( pkg->GetPSideCluster(0)->GetTotalSignal(),
1736 pkg->GetNSideCluster(0)->GetTotalSignal() );
e8189707 1737 for (Int_t i = 0; i< pkg->GetNumOfClustersP(); i++)
b0f5e3fc 1738 {
1739
1740 p = pkg->GetPSideClusterIdx(i);
1741 curPcl= GetPSideCluster(p);
1742 nc=curPcl->GetCrossNo();
1743 ps=curPcl->GetTotalSignal();
1744
e8189707 1745 for (Int_t j = 0; j< nc; j++)
b0f5e3fc 1746 {
1747 n=curPcl->GetCross(j);
1748 curNcl= GetNSideCluster(n);
1749 ns=curNcl->GetTotalSignal();
1750 chi = DistToPML(ps, ns);
1751
1752 if (chi <= bestchi)
1753 {
1754 bestchi = chi;
1755 in = n;
1756 ip = p;
1757 }
1758 }
1759 }
1760
1761}
1762
1763
1764//------------------------------------------------------
1765void AliITSClusterFinderSSD::CalcStepFactor(Float_t Psteo, Float_t Nsteo )
1766{
1767
1768
1769//Piotr Krzysztof Skowronski
1770//Warsaw University of Technology
1771//skowron@if.pw.edu.pl
1772
1773
1774 // 95 is the pitch, 4000 - dimension along z ?
1775 /*
1776 fSFF = ( (Int_t) (Psteo*40000/95 ) );// +1;
1777 fSFB = ( (Int_t) (Nsteo*40000/95 ) );// +1;
1778 */
1779
1780
1781 Float_t dz=fSegmentation->Dz();
1782
1783 fSFF = ( (Int_t) (Psteo*dz/fPitch ) );// +1;
1784 fSFB = ( (Int_t) (Nsteo*dz/fPitch ) );// +1;
1785
1786}
1787
1788
1789//-----------------------------------------------------------
1790AliITSclusterSSD* AliITSClusterFinderSSD::GetPSideCluster(Int_t idx)
1791{
1792
1793
1794//Piotr Krzysztof Skowronski
1795//Warsaw University of Technology
1796//skowron@if.pw.edu.pl
1797
1798
1799
1800 if((idx<0)||(idx>=fNClusterP))
1801 {
1802 printf("AliITSClusterFinderSSD::GetPSideCluster : index out of range\n");
1803 return 0;
1804 }
1805 else
1806 {
1807 return (AliITSclusterSSD*)((*fClusterP)[idx]);
1808 }
1809}
1810
1811//-------------------------------------------------------
1812AliITSclusterSSD* AliITSClusterFinderSSD::GetNSideCluster(Int_t idx)
1813{
1814
1815
1816//Piotr Krzysztof Skowronski
1817//Warsaw University of Technology
1818//skowron@if.pw.edu.pl
1819
1820 if((idx<0)||(idx>=fNClusterN))
1821 {
1822 printf("AliITSClusterFinderSSD::GetNSideCluster : index out of range\n");
1823 return 0;
1824 }
1825 else
1826 {
1827 return (AliITSclusterSSD*)((*fClusterN)[idx]);
1828 }
1829}
1830
1831//--------------------------------------------------------
1832AliITSclusterSSD* AliITSClusterFinderSSD::GetCluster(Int_t idx, Bool_t side)
1833{
1834
1835
1836//Piotr Krzysztof Skowronski
1837//Warsaw University of Technology
1838//skowron@if.pw.edu.pl
1839
1840 return (side) ? GetPSideCluster(idx) : GetNSideCluster(idx);
1841}
1842
1843
1844//--------------------------------------------------------
1845void AliITSClusterFinderSSD::ConsumeClusters()
1846{
1847
e8189707 1848
1849 for ( Int_t i=0;i<fNPackages;i++)
b0f5e3fc 1850 {
1851 ((AliITSpackageSSD*)((*fPackages)[i]))->ConsumeClusters();
1852 }
1853
1854}
1855
1856
1857//--------------------------------------------------------
1858void AliITSClusterFinderSSD::ReconstructNotConsumedClusters()
1859{
1860 Int_t i;
1861 AliITSclusterSSD *cluster;
1862 Float_t pos;
1863 Float_t sig;
1864 Float_t sigerr;
1865 Float_t x1,x2,z1,z2;
1866
1867 Float_t Dx = fSegmentation->Dx();
1868 Float_t Dz = fSegmentation->Dz();
1869
1870 const Float_t kdEdXtoQ = 2.778e+8; // GeV -> number of e-hole pairs
1871 const Float_t kconv = 1.0e-4;
1872 const Float_t kRMSx = 20.0*kconv; // microns->cm ITS TDR Table 1.3
1873 const Float_t kRMSz = 830.0*kconv; // microns->cm ITS TDR Table 1.3
1874
e8189707 1875 Int_t stripP,stripN;
1876 AliITSdigitSSD *dig;
1877 for (i=0;i<fNClusterP;i++)
b0f5e3fc 1878 {
b0f5e3fc 1879 cluster = GetPSideCluster(i);
1880 if (!cluster->IsConsumed())
1881 {
1882 if( (pos = cluster->GetPosition()) < fSFB)
1883 {
1884 sig = cluster->GetTotalSignal();
1885
be33dccb 1886 sig += sig/fPNsignalRatio;
b0f5e3fc 1887
1888 sigerr = cluster->GetTotalSignalError();
1889 x1 = -Dx/2 + pos *fPitch;
1890 z1 = -Dz/2;
1891
1892 x2 = pos;
1893 z2 = 1;
1894 GetCrossing (x2,z2);
1895 AliITSRawClusterSSD cnew;
1896 Int_t nstripsP=cluster->GetNumOfDigits();
1897 cnew.fMultiplicity=nstripsP;
1898 cnew.fMultiplicityN=0;
1899 /*
1900 if (nstripsP > 100) {
1901 printf("multiplicity > 100 - increase the dimension of the arrays %d\n",nstripsP);
1902 nstripsP=100;
1903 }
e8189707 1904 for(int k=0;k<nstripsP;k++) {
b0f5e3fc 1905 // check if 'clusterP->GetDigitStripNo(i)' returns
1906 // the digit index and not smth else
1907 cnew.fIndexMap[k] = cluster->GetDigitStripNo(k);
1908 }
1909 */
1910 cnew.fQErr=sigerr;
1911 //cnew.fProbability=0.75;
1912 fITS->AddCluster(2,&cnew);
e8189707 1913 fSegmentation->GetPadIxz((x1+x2)/2,(z1+z2)/2,stripP,stripN);
1914 dig = (AliITSdigitSSD*)fMap->GetHit(1,stripP);
b0f5e3fc 1915 // add the rec point info
1916 AliITSRecPoint rnew;
1917 rnew.SetX(kconv*(x1+x2)/2);
1918 rnew.SetZ(kconv*(z1+z2)/2);
1919 rnew.SetQ(sig);
1920 rnew.SetdEdX(sig/kdEdXtoQ);
1921 rnew.SetSigmaX2( kRMSx* kRMSx);
1922 rnew.SetSigmaZ2( kRMSz* kRMSz);
1923 rnew.SetProbability(0.75);
e8189707 1924 rnew.fTracks[0]=dig->fTracks[0];
1925 rnew.fTracks[1]=dig->fTracks[1];
1926 rnew.fTracks[2]=dig->fTracks[2];
b0f5e3fc 1927 fITS->AddRecPoint(rnew);
1928 /*
1929 fPointsM->AddLast( (TObject*)
1930 new AliITSRecPointSSD((x1+x2)/2 ,0,0,0,(z1+z2)/2,0,sig,sigerr,fLayer,fLad,fDet,cluster, 0.75) );
1931 */
1932
1933 if(debug) cout<<"\n"<<": SINGLE SIDE ImpactPoint Created: X = "
1934 <<(x1+x2)/2<<" Z = "<<(z1+z2)/2<<"; Pos = "<<pos;
1935 }
1936 }
1937 }
1938
e8189707 1939 for (i=0;i<fNClusterN;i++)
b0f5e3fc 1940 {
b0f5e3fc 1941 cluster = GetNSideCluster(i);
1942 if (!cluster->IsConsumed())
1943 {
1944 if( (pos = cluster->GetPosition()) < fSFF)
1945 {
1946 sig = cluster->GetTotalSignal();
1947
be33dccb 1948 sig += sig*fPNsignalRatio;
b0f5e3fc 1949
1950 sigerr = cluster->GetTotalSignalError();
1951
1952 x1 = -Dx/2 + pos *fPitch;
1953 z1 = Dz/2;
1954
1955 x2 = 1;
1956 z2 = pos;
1957
1958 GetCrossing (x2,z2);
1959 AliITSRawClusterSSD cnew;
1960 Int_t nstripsN=cluster->GetNumOfDigits();
1961 cnew.fMultiplicity=0;
1962 cnew.fMultiplicityN=nstripsN;
1963 /*
1964 if (nstripsN > 100) {
1965 printf("multiplicity > 100 - increase the dimension of the arrays %d\n",nstripsN);
1966 nstripsN=100;
1967 }
e8189707 1968 for(int k=0;k<nstripsN;k++) {
b0f5e3fc 1969 // check if 'clusterP->GetDigitStripNo(i)' returns
1970 // the digit index and not smth else
1971 cnew.fIndexMapN[k] = cluster->GetDigitStripNo(k);
1972 }
1973 */
1974 cnew.fQErr=sigerr;
1975 //cnew.fProbability=0.75;
1976 fITS->AddCluster(2,&cnew);
1977 // add the rec point info
e8189707 1978 fSegmentation->GetPadIxz((x1+x2)/2,(z1+z2)/2,stripP,stripN);
1979 dig = (AliITSdigitSSD*)fMap->GetHit(0,stripN);
b0f5e3fc 1980 AliITSRecPoint rnew;
1981 rnew.SetX(kconv*(x1+x2)/2);
1982 rnew.SetZ(kconv*(z1+z2)/2);
1983 rnew.SetQ(sig);
1984 rnew.SetdEdX(sig/kdEdXtoQ);
1985 rnew.SetSigmaX2( kRMSx* kRMSx);
1986 rnew.SetSigmaZ2( kRMSz* kRMSz);
1987 rnew.SetProbability(0.75);
e8189707 1988 rnew.fTracks[0]=dig->fTracks[0];
1989 rnew.fTracks[1]=dig->fTracks[1];
1990 rnew.fTracks[2]=dig->fTracks[2];
b0f5e3fc 1991 fITS->AddRecPoint(rnew);
1992 /*
1993 fPointsM->AddLast( (TObject*)
1994 new AliITSRecPointSSD((x1+x2)/2 ,0,0,0,(z1+z2)/2,0,sig,sigerr,fLayer,fLad,fDet,cluster, 0.75) );
1995 */
1996
1997 if(debug) cout<<"\n"<<": SINGLE SIDE ImpactPoint Created: X = "
1998 <<(x1+x2)/2<<" Z = "<<(z1+z2)/2<<"; Pos = "<<pos;
1999
2000 }
2001 }
2002 }
2003
2004
2005}
2006
2007//_______________________________________________________________________
2008
2009Bool_t AliITSClusterFinderSSD::GetCrossing (Float_t &P, Float_t &N)
2010{
2011
2012 Float_t Dx = fSegmentation->Dx();
2013 Float_t Dz = fSegmentation->Dz();
2014
2015 //P==X
2016 //N==Z
2017
2018 Float_t dx = 0.1;
2019 P *= fPitch;
2020 N *= fPitch;
2021
2022 //P = (kZ * fTan + N + P)/2.0; // x coordinate
2023 //N = kZ - (P-N)/fTan; // z coordinate
2024
2025 if(fTanP + fTanN == 0) return kFALSE;
2026
2027 N = (N - P + fTanN * Dz) / (fTanP + fTanN); // X coordinate
2028 P = P + fTanP * N; // Y coordinate
2029
2030 P -= Dx/2;
2031 N -= Dz/2;
2032
2033 // N = -(N - P + kZ/2*(fTanP + fTanN))/(fTanP + fTanN);
2034 // P = (fTanP*(N-kX/2-kZ*fTanN/2) + fTanN*(P-kX/2-kZ*fTanP/2) )/(fTanP + fTanN);
2035
2036 if ( ( N < -(Dz/2+dx) ) || ( N > (Dz/2+dx) ) ) return kFALSE;
2037 if ( ( P < -(Dx/2+dx) ) || ( P > (Dx/2+dx) ) ) return kFALSE;
2038
2039 return kTRUE;
2040}
2041
2042
2043//_________________________________________________________________________
2044
2045void AliITSClusterFinderSSD::GetCrossingError(Float_t& dP, Float_t& dN)
2046{
2047 Float_t dz, dx;
2048
2049 dz = TMath::Abs(( dP + dN )*fPitch/(fTanP + fTanN) );
2050 dx = fPitch*(TMath::Abs(dP*(1 - fTanP/(fTanP + fTanN))) +
2051 TMath::Abs(dN *fTanP/(fTanP + fTanN) ));
2052
2053 dN = dz;
2054 dP = dx;
2055}
2056
2057