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