2351d7362df69be9c10d5e97fd196c55e9086117
[u/mrichter/AliRoot.git] / MUON / AliMUONSegmentationV01.cxx
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 $Log$
18 Revision 1.4  2000/07/03 11:54:57  morsch
19 AliMUONSegmentation and AliMUONHitMap have been replaced by AliSegmentation and AliHitMap in STEER
20 The methods GetPadIxy and GetPadXxy of AliMUONSegmentation have changed name to GetPadI and GetPadC.
21
22 Revision 1.3  2000/06/29 12:34:09  morsch
23 AliMUONSegmentation class has been made independent of AliMUONChamber. This makes
24 it usable with any other geometry class. The link to the object to which it belongs is
25 established via an index. This assumes that there exists a global geometry manager
26 from which the pointer to the parent object can be obtained (in our case gAlice).
27
28 Revision 1.2  2000/06/15 07:58:48  morsch
29 Code from MUON-dev joined
30
31 Revision 1.1.2.1  2000/06/09 21:37:30  morsch
32 AliMUONSegmentationV01 code  from  AliMUONSegResV01.cxx
33
34 */
35
36
37 /////////////////////////////////////////////////////
38 //  Segmentation and Response classes version 01   //
39 /////////////////////////////////////////////////////
40
41 #include <TBox.h> 
42 #include <TF1.h> 
43 #include <TObjArray.h>
44 #include <iostream.h>
45
46 #include "AliMUONSegmentationV01.h"
47 #include "AliMUON.h"
48
49
50
51 //___________________________________________
52 ClassImp(AliMUONSegmentationV01)
53
54 AliMUONSegmentationV01::AliMUONSegmentationV01(const AliMUONSegmentationV01& segmentation)
55 {
56 // Dummy copy constructor
57 }
58 AliMUONSegmentationV01::AliMUONSegmentationV01() 
59 {
60 // Default constructor
61     fNsec=4;
62     fRSec.Set(fNsec);    
63     fNDiv.Set(fNsec);      
64     fDpxD.Set(fNsec);      
65     fRSec[0]=fRSec[1]=fRSec[2]=fRSec[3]=0;     
66     fNDiv[0]=fNDiv[1]=fNDiv[2]=fNDiv[3]=0;     
67     fDpxD[0]=fDpxD[1]=fDpxD[2]=fDpxD[3]=0;     
68     fCorr = new TObjArray(3);
69     (*fCorr)[0]=0;
70     (*fCorr)[1]=0;
71     (*fCorr)[2]=0;
72
73
74 Float_t AliMUONSegmentationV01::Dpx(Int_t isec)
75 {
76 //
77 // Returns x-pad size for given sector isec
78    return fDpxD[isec];
79 }
80
81 Float_t AliMUONSegmentationV01::Dpy(Int_t isec)
82 {
83 //
84 // Returns y-pad size for given sector isec
85    return fDpy;
86 }
87     
88 void   AliMUONSegmentationV01::SetSegRadii(Float_t  r[4])
89 {
90 //
91 // Set the radii of the segmentation zones 
92     for (Int_t i=0; i<4; i++) {
93         fRSec[i]=r[i];
94         printf("\n R %d %f \n",i,fRSec[i]);
95         
96     }
97 }
98
99
100 void AliMUONSegmentationV01::SetPadDivision(Int_t ndiv[4])
101 {
102 //
103 // Defines the pad size perp. to the anode wire (y) for different sectors. 
104 // Pad sizes are defined as integral fractions ndiv of a basis pad size
105 // fDpx
106 // 
107     for (Int_t i=0; i<4; i++) {
108         fNDiv[i]=ndiv[i];
109         printf("\n Ndiv %d %d \n",i,fNDiv[i]);
110     }
111     ndiv[0]=ndiv[1];
112 }
113
114
115 void AliMUONSegmentationV01::Init(Int_t chamber)
116 {
117 //
118 //  Fill the arrays fCx (x-contour) and fNpxS (ix-contour) for each sector
119 //  These arrays help in converting from real to pad co-ordinates and
120 //  vice versa.
121 //  This version approximates concentric segmentation zones
122 //
123     Int_t isec;
124     printf("\n Initialise segmentation v01 -- test !!!!!!!!!!!!!! \n");
125     fNpy=Int_t(fRSec[fNsec-1]/fDpy)+1;
126
127     fDpxD[fNsec-1]=fDpx;
128     if (fNsec > 1) {
129         for (Int_t i=fNsec-2; i>=0; i--){
130             fDpxD[i]=fDpxD[fNsec-1]/fNDiv[i];
131             printf("\n test ---dx %d %f \n",i,fDpxD[i]);
132         }
133     }
134 //
135 // fill the arrays defining the pad segmentation boundaries
136     Float_t ry;
137     Int_t   dnx;
138     Int_t   add;
139 //
140 //  loop over sections
141     for(isec=0; isec<fNsec; isec++) {
142 //  
143 //  loop over pads along the aode wires
144         for (Int_t iy=1; iy<=fNpy; iy++) {
145 //
146             Float_t x=iy*fDpy-fDpy/2;
147             if (x > fRSec[isec]) {
148                 fNpxS[isec][iy]=0;
149                 fCx[isec][iy]=0;
150             } else {
151                 ry=TMath::Sqrt(fRSec[isec]*fRSec[isec]-x*x);
152                 if (isec > 1) {
153                     dnx= Int_t((ry-fCx[isec-1][iy])/fDpxD[isec]);
154                     if (isec < fNsec-1) {
155                         if (TMath::Odd((Long_t)dnx)) dnx++;             
156                     }
157                     fNpxS[isec][iy]=fNpxS[isec-1][iy]+dnx;
158                     fCx[isec][iy]=fCx[isec-1][iy]+dnx*fDpxD[isec];
159                 } else if (isec == 1) {
160                     dnx= Int_t((ry-fCx[isec-1][iy])/fDpxD[isec]);
161                     fNpxS[isec][iy]=fNpxS[isec-1][iy]+dnx;
162                     add=4 - (fNpxS[isec][iy])%4;
163                     if (add < 4) fNpxS[isec][iy]+=add; 
164                     dnx=fNpxS[isec][iy]-fNpxS[isec-1][iy];
165                     fCx[isec][iy]=fCx[isec-1][iy]+dnx*fDpxD[isec];
166                 } else {
167                     dnx=Int_t(ry/fDpxD[isec]);
168                     fNpxS[isec][iy]=dnx;
169                     fCx[isec][iy]=dnx*fDpxD[isec];
170                 }
171             }
172         } // y-pad loop
173     } // sector loop
174 }
175
176 Int_t AliMUONSegmentationV01::Sector(Int_t ix, Int_t iy)
177 {
178 // Returns sector number for given pad position
179 //
180     Int_t absix=TMath::Abs(ix);
181     Int_t absiy=TMath::Abs(iy);
182     Int_t isec=0;
183     for (Int_t i=0; i<fNsec; i++) {
184         if (absix<=fNpxS[i][absiy]){
185             isec=i;
186             break;
187         }
188     }
189     return isec;
190 }
191
192 void AliMUONSegmentationV01::
193 GetPadI(Float_t x, Float_t y, Int_t &ix, Int_t &iy)
194 {
195 //  Returns pad coordinates (ix,iy) for given real coordinates (x,y)
196 //
197     iy = (y>0)? Int_t(y/fDpy)+1 : Int_t(y/fDpy)-1;
198     if (iy >  fNpy) iy= fNpy;
199     if (iy < -fNpy) iy=-fNpy;
200 //
201 //  Find sector isec
202     Int_t isec=-1;
203     Float_t absx=TMath::Abs(x);
204     Int_t absiy=TMath::Abs(iy);
205     for (Int_t i=0; i < fNsec; i++) {
206         if (absx <= fCx[i][absiy]) {
207             isec=i;
208             break;
209         }
210     }
211     if (isec>0) {
212         ix= Int_t((absx-fCx[isec-1][absiy])/fDpxD[isec])
213             +fNpxS[isec-1][absiy]+1;
214     } else if (isec == 0) {
215         ix= Int_t(absx/fDpxD[isec])+1;
216     } else {
217         ix=fNpxS[fNsec-1][absiy]+1;     
218     }
219     ix = (x>0) ? ix:-ix;
220 }
221
222 void AliMUONSegmentationV01::
223 GetPadC(Int_t ix, Int_t iy, Float_t &x, Float_t &y)
224 {
225 //  Returns real coordinates (x,y) for given pad coordinates (ix,iy)
226 //
227     y = (iy>0) ? Float_t(iy*fDpy)-fDpy/2. : Float_t(iy*fDpy)+fDpy/2.;
228 //
229 //  Find sector isec
230     Int_t isec=AliMUONSegmentationV01::Sector(ix,iy);
231 //
232     Int_t absix=TMath::Abs(ix);
233     Int_t absiy=TMath::Abs(iy);
234     if (isec) {
235         x=fCx[isec-1][absiy]+(absix-fNpxS[isec-1][absiy])*fDpxD[isec];
236         x=(ix>0) ?  x-fDpxD[isec]/2 : -x+fDpxD[isec]/2;
237     } else {
238         x=y=0;
239     }
240 }
241
242 void AliMUONSegmentationV01::
243 SetPad(Int_t ix, Int_t iy)
244 {
245     //
246     // Sets virtual pad coordinates, needed for evaluating pad response 
247     // outside the tracking program 
248     GetPadC(ix,iy,fx,fy);
249     fSector=Sector(ix,iy);
250 }
251
252
253 void AliMUONSegmentationV01::
254 FirstPad(Float_t xhit, Float_t yhit, Float_t dx, Float_t dy)
255 {
256 // Initialises iteration over pads for charge distribution algorithm
257 //
258     //
259     // Find the wire position (center of charge distribution)
260     Float_t x0a=GetAnod(xhit);
261     fxhit=x0a;
262     fyhit=yhit;
263     
264     //
265     // and take fNsigma*sigma around this center
266     Float_t x01=x0a  - dx;
267     Float_t x02=x0a  + dx;
268     Float_t y01=yhit - dy;
269     Float_t y02=yhit + dy;
270     //
271     // find the pads over which the charge distributes
272     GetPadI(x01,y01,fixmin,fiymin);
273     GetPadI(x02,y02,fixmax,fiymax);
274     fxmin=x01;
275     fxmax=x02;
276     fymin=y01;
277     fymax=y02;
278     
279     // 
280     // Set current pad to lower left corner
281     if (fixmax < fixmin) fixmax=fixmin;
282     if (fiymax < fiymin) fiymax=fiymin;    
283     fix=fixmin;
284     fiy=fiymin;
285     GetPadC(fix,fiy,fx,fy);
286 }
287
288
289 void AliMUONSegmentationV01::NextPad()
290 {
291 // Stepper for the iteration over pads
292 //
293 // Step to next pad in the integration region
294   // 
295   // Step to next pad in integration region
296     Float_t xc,yc;
297     Int_t   iyc;
298     
299 //  step from left to right    
300     if (fx < fxmax && fx != 0) {
301         if (fix==-1) fix++;
302         fix++;
303 //  step up 
304     } else if (fiy != fiymax) {
305         if (fiy==-1) fiy++;
306         fiy++;
307 //      get y-position of next row (yc), xc not used here       
308         GetPadC(fix,fiy,xc,yc);
309 //      get x-pad coordiante for first pad in row (fix)
310         GetPadI(fxmin,yc,fix,iyc);
311     } else {
312         printf("\n Error: Stepping outside integration region\n ");
313     }
314     GetPadC(fix,fiy,fx,fy);
315     fSector=Sector(fix,fiy);
316     if (MorePads() && 
317         (fSector ==-1 || fSector==0)) 
318         NextPad();
319 }
320
321 Int_t AliMUONSegmentationV01::MorePads()
322 // Stopping condition for the iterator over pads
323 //
324 // Are there more pads in the integration region
325 {
326     if ((fx >= fxmax  && fiy >= fiymax) || fy==0) {
327         return 0;
328     } else {
329         return 1;
330     }
331 }
332
333 void AliMUONSegmentationV01::
334 IntegrationLimits(Float_t& x1,Float_t& x2,Float_t& y1, Float_t& y2)
335 {
336 //  Returns integration limits for current pad
337 //
338     x1=fxhit-fx-Dpx(fSector)/2.;
339     x2=x1+Dpx(fSector);
340     y1=fyhit-fy-Dpy(fSector)/2.;
341     y2=y1+Dpy(fSector);    
342 }
343
344 void AliMUONSegmentationV01::
345 Neighbours(Int_t iX, Int_t iY, Int_t* Nlist, Int_t Xlist[10], Int_t Ylist[10])
346 {
347 // Returns list of next neighbours for given Pad (iX, iY)
348 //
349     const Float_t kEpsilon=fDpy/1000;
350     
351     Float_t x,y;
352     Int_t   ixx, iyy, isec1;
353 //
354     Int_t   isec0=AliMUONSegmentationV01::Sector(iX,iY);
355     Int_t i=0;
356 //    
357 //  step right
358     Xlist[i]=iX+1;
359     if (Xlist[i]==0) Xlist[i]++;
360     Ylist[i++]=iY;
361 //
362 //  step left    
363     Xlist[i]=iX-1;
364     if (Xlist[i]==0) Xlist[i]--;
365     Ylist[i++]=iY;
366 //
367 //  step up 
368     AliMUONSegmentationV01::GetPadC(iX,iY,x,y);
369     AliMUONSegmentationV01::GetPadI(x+kEpsilon,y+fDpy,ixx,iyy);
370     Xlist[i]=ixx;
371     Ylist[i++]=iyy;
372     isec1=AliMUONSegmentationV01::Sector(ixx,iyy);
373     if (isec1==isec0) {
374 //
375 //  no sector boundary crossing
376 //      Xlist[i]=ixx+1;
377 //      Ylist[i++]=iY+1;
378         
379 //      Xlist[i]=ixx-1;
380 //      Ylist[i++]=iY+1;
381     } else if (isec1 < isec0) {
382 // finer segmentation
383 //      Xlist[i]=ixx+1;
384 //      Ylist[i++]=iY+1;
385         
386         Xlist[i]=ixx-1;
387         Ylist[i++]=iyy;
388         
389 //      Xlist[i]=ixx-2;
390 //      Ylist[i++]=iY+1;
391     } else {
392 // coarser segmenation  
393 /*
394         if (TMath::Odd(iX-fNpxS[isec1-1][iY+1])) {
395             Xlist[i]=ixx-1;
396             Ylist[i++]=iY+1;
397         } else {
398             Xlist[i]=ixx+1;
399             Ylist[i++]=iY+1;
400         }
401 */
402     }
403
404 //
405 // step down 
406     AliMUONSegmentationV01::GetPadC(iX,iY,x,y);
407     AliMUONSegmentationV01::GetPadI(x+kEpsilon,y-fDpy,ixx,iyy);
408     Xlist[i]=ixx;
409     Ylist[i++]=iyy;
410     isec1=AliMUONSegmentationV01::Sector(ixx,iyy);
411     if (isec1==isec0) {
412 //
413 //  no sector boundary crossing
414 /*
415     Xlist[i]=ixx+1;
416     Ylist[i++]=iY-1;
417         
418     Xlist[i]=ixx-1;
419     Ylist[i++]=iY-1;
420 */
421     } else if (isec1 < isec0) {
422 // finer segmentation
423 //      Xlist[i]=ixx+1;
424 //      Ylist[i++]=iY-1;
425         
426         Xlist[i]=ixx-1;
427         Ylist[i++]=iyy;
428         
429 //      Xlist[i]=ixx-2;
430 //      Ylist[i++]=iY-1;
431     } else {
432 // coarser segmentation 
433 /*
434         if (TMath::Odd(iX-fNpxS[isec1-1][iY-1])) {
435             Xlist[i]=ixx-1;
436             Ylist[i++]=iY-1;
437         } else {
438             Xlist[i]=ixx+1;
439             Ylist[i++]=iY-1;
440         }
441 */
442     }
443     *Nlist=i;
444 }
445
446 void AliMUONSegmentationV01::GiveTestPoints(Int_t &n, Float_t *x, Float_t *y)
447 {
448 // Returns test point on the pad plane.
449 // Used during determination of the segmoid correction of the COG-method
450
451     n=3;
452     x[0]=(fRSec[0]+fRSec[1])/2/TMath::Sqrt(2.);
453     y[0]=x[0];
454     x[1]=(fRSec[1]+fRSec[2])/2/TMath::Sqrt(2.);
455     y[1]=x[1];
456     x[2]=(fRSec[2]+fRSec[3])/2/TMath::Sqrt(2.);
457     y[2]=x[2];
458 }
459
460 void AliMUONSegmentationV01::Draw(const char *)
461 {
462 // Draws the segmentation zones
463 //
464     TBox *box;
465     
466     Float_t dx=0.95/fCx[3][1]/2;
467     Float_t dy=0.95/(Float_t(Npy()))/2;
468     Float_t x0,y0,x1,y1;
469     Float_t xc=0.5;
470     Float_t yc=0.5;
471     
472     for (Int_t iy=1; iy<Npy(); iy++)
473     {
474         for (Int_t isec=0; isec<4; isec++) {
475             if (isec==0) {
476                 x0=0;
477                 x1=fCx[isec][iy]*dx;
478             } else {
479                 x0=fCx[isec-1][iy]*dx;
480                 x1=fCx[isec][iy]*dx;
481             }
482             y0=Float_t(iy-1)*dy;
483             y1=y0+dy;
484             box=new TBox(x0+xc,y0+yc,x1+xc,y1+yc);
485             box->SetFillColor(isec+1);
486             box->Draw();
487
488             box=new TBox(-x1+xc,y0+yc,-x0+xc,y1+yc);
489             box->SetFillColor(isec+1);
490             box->Draw();
491
492             box=new TBox(x0+xc,-y1+yc,x1+xc,-y0+yc);
493             box->SetFillColor(isec+1);
494             box->Draw();
495
496             box=new TBox(-x1+xc,-y1+yc,-x0+xc,-y0+yc);
497             box->SetFillColor(isec+1);
498             box->Draw();
499         }
500     }
501 }
502 void AliMUONSegmentationV01::SetCorrFunc(Int_t isec, TF1* func)
503 {
504     (*fCorr)[isec]=func;
505 }
506
507 TF1* AliMUONSegmentationV01::CorrFunc(Int_t isec)
508
509     return (TF1*) (*fCorr)[isec];
510 }
511
512 AliMUONSegmentationV01& AliMUONSegmentationV01::operator 
513 =(const AliMUONSegmentationV01 & rhs)
514 {
515 // Dummy assignment operator
516     return *this;
517 }