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