8e1c6b28ae97697a39a0c3a1f6efe7dd6eb1c247
[u/mrichter/AliRoot.git] / MUON / AliMUONSegmentationV1.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.3  2000/06/29 12:34:09  morsch
19 AliMUONSegmentation class has been made independent of AliMUONChamber. This makes
20 it usable with any other geometry class. The link to the object to which it belongs is
21 established via an index. This assumes that there exists a global geometry manager
22 from which the pointer to the parent object can be obtained (in our case gAlice).
23
24 Revision 1.2  2000/06/15 07:58:48  morsch
25 Code from MUON-dev joined
26
27 Revision 1.1.2.2  2000/06/12 07:57:43  morsch
28 include TMath.cxx
29
30 Revision 1.1.2.1  2000/06/09 21:41:29  morsch
31 AliMUONSegmentationV1 code  from  AliMUONSegResV1.cxx
32
33 */
34
35
36 /////////////////////////////////////////////////////////
37 //  Manager and hits classes for set:MUON version LYON //
38 /////////////////////////////////////////////////////////
39
40 #include <TMath.h>
41 #include "AliMUONChamber.h"
42 #include "AliMUONSegmentationV1.h"
43 #include "AliRun.h"
44 #include "AliMUON.h"
45
46 //___________________________________________
47 ClassImp(AliMUONSegmentationV1)
48
49 AliMUONSegmentationV1::AliMUONSegmentationV1(const AliMUONSegmentationV1& segmentation)
50 {
51 // Dummy copy constructor
52 }
53
54
55 AliMUONSegmentationV1::AliMUONSegmentationV1()
56
57 {
58   // initizalize the class with default settings
59     fNzone=1;
60     fDAnod=0.0;
61     fDpx=0.0;
62     fDpx=0.0; // forces crash if not initialized by user
63     fNZoneCut[0]=0;
64     fSensOffset=0;
65 }
66
67
68 void AliMUONSegmentationV1::Init(Int_t chamber)
69 {
70     // valid only for T5/6
71     // beware : frMin is SENSITIVE radius by definition.
72     AliMUON *pMUON  = (AliMUON *) gAlice->GetModule("MUON");
73     AliMUONChamber* iChamber=&(pMUON->Chamber(chamber));
74
75     frSensMin2 = (iChamber->RInner())*(iChamber->RInner());
76     frSensMax2 = (iChamber->ROuter())*(iChamber->ROuter());
77     fNpx=(Int_t) (iChamber->ROuter()/fDpx) + 1;
78     fNpy=(Int_t) (iChamber->ROuter()/fDpy) + 1;
79     //    fNwire=3;
80     DefaultCut();
81     fCorr=0;
82 }
83
84 void AliMUONSegmentationV1::DefaultCut(void)
85 {
86 // Set the default cuts
87     SetNzone(3);
88     AddCut(0,5*6,18*8);
89     AddCut(0,9*6,15*8);
90     AddCut(0,11*6,12*8);
91     AddCut(0,12*6,9*8);
92     AddCut(0,13*6,6*8);
93     AddCut(1,6*6,20*12);
94     AddCut(1,12*6,18*12);
95     AddCut(1,15*6,15*12);
96     AddCut(1,18*6,12*12);
97     AddCut(1,21*6,9*12);
98     SetSensOffset(3.0);
99     SetDAnod(0.325);
100 }
101
102 Int_t AliMUONSegmentationV1::GetiAnod(Float_t xhit)
103 {
104 // Get anode number
105     Int_t kwire=Int_t((TMath::Abs(xhit)-fSensOffset)/fDAnod)+1;
106     return (xhit>0) ? kwire : -kwire ;
107 }
108
109 Float_t AliMUONSegmentationV1::GetAnod(Float_t xhit)
110 {
111 // Get anode position
112     Int_t kwire=Int_t((TMath::Abs(xhit)-fSensOffset)/fDAnod)+1; // to be compatible ...
113     return (xhit>0) ? fDAnod*(kwire-0.5)+fSensOffset : -fDAnod*(kwire-0.5)-fSensOffset ;
114 }
115
116
117 void AliMUONSegmentationV1::SetPadSize(Float_t p1, Float_t p2)
118 {
119 // For chamber T5/6 p1 and p2 should be same for each zone
120     fDpx=p1;
121     fDpy=p2;
122 }
123
124 void AliMUONSegmentationV1::
125 GetPadI(Float_t x, Float_t y, Int_t &ix, Int_t &iy)
126 {
127 //  returns pad coordinates (ix,iy) for given real coordinates (x,y)
128 //
129     ix = (x>0)? Int_t((x-fSensOffset)/fDpx)+1 : Int_t((x+fSensOffset)/fDpx)-1;
130     iy = (y>0)? Int_t((y-fSensOffset)/fDpy)+1 : Int_t((y+fSensOffset)/fDpy)-1;
131 }
132
133 void AliMUONSegmentationV1::
134 GetPadC(Int_t ix, Int_t iy, Float_t &x, Float_t &y)
135 {
136 //  returns real coordinates (x,y) for given pad coordinates (ix,iy)
137 //
138     x = (ix>0) ? (Float_t(ix)-0.5)*fDpx+fSensOffset : (Float_t(ix)+0.5)*fDpx-fSensOffset;
139     y = (iy>0) ? (Float_t(iy)-0.5)*fDpy+fSensOffset : (Float_t(iy)+0.5)*fDpy-fSensOffset;
140 }
141
142 void AliMUONSegmentationV1::AddCut(Int_t Zone, Int_t nX, Int_t nY)
143 {
144 // the pad nX,nY is last INSIDE zone Zone. First pad is labelled 1 and not 0
145     if (Zone+1>=fNzone) 
146 // no cut for last Zone : it is the natural boundary of the chamber
147         printf("AliMUONSegmentationV1::AddCut ==> Zone %d not allowed !\n",Zone);
148     fZoneX[Zone][fNZoneCut[Zone]] = nX;
149     fZoneY[Zone][fNZoneCut[Zone]] = nY;
150     fNZoneCut[Zone]++;
151 }
152
153 Int_t AliMUONSegmentationV1::GetZone(Float_t X, Float_t Y)
154 {
155 // Get segmentation zone
156     Int_t iX, iY;
157     GetPadI(X,Y,iX,iY);
158     return GetZone( iX , iY );
159 }
160
161 Int_t AliMUONSegmentationV1::GetZone(Int_t nX, Int_t nY)
162 {
163 // Beware : first pad begins at 1 !!
164     Int_t aX =  TMath::Abs(nX);
165     Int_t aY =  TMath::Abs(nY);
166     Int_t zone=fNzone-1;
167     for (Int_t iZone=fNzone-2;iZone>=0;iZone--) 
168     {
169         for (Int_t iCut=0;iCut<fNZoneCut[iZone];iCut++)
170             if ( aY<=fZoneY[iZone][iCut] && aX<=fZoneX[iZone][iCut] )
171             {
172                 zone=iZone;
173                 break;
174             } 
175     }
176     return zone;
177 }
178
179 void AliMUONSegmentationV1::
180 SetHit(Float_t xhit, Float_t yhit)
181 {
182 // Find the wire position (center of charge distribution)
183     fxhit=xhit;
184     fyhit=yhit;
185 }
186
187 void AliMUONSegmentationV1::
188 SetPad(Int_t ix, Int_t iy)
189 {
190 // Set current pad position
191     GetPadC(ix,iy,fx,fy);
192 }
193
194
195 void AliMUONSegmentationV1::SetPadCoord(Int_t iX, Int_t iY)
196 {    
197 // Set current pad coordinates
198 GetPadC(iX,iY,fx,fy);
199  Float_t radius2;
200  if ( ( (radius2=fx*fx+fy*fy) > frSensMax2 || radius2 < frSensMin2 ) 
201       && MorePads() )
202      NextPad();
203 }
204
205 void AliMUONSegmentationV1::FirstPad(Float_t xhit, Float_t yhit, Float_t dx, Float_t dy)
206 {
207     //
208     // Find the wire position (center of charge distribution)
209     Float_t x0a=GetAnod(xhit);
210     fxhit=x0a;
211     fyhit=yhit;
212     //
213     // and take fNsigma*sigma around this center
214     Float_t x01=x0a  - dx;
215     Float_t x02=x0a  + dx;
216     Float_t y01=yhit - dy;
217     Float_t y02=yhit + dy;
218
219     // Do not cross over frames...
220     if (x01 * x0a < 0) 
221       x01 = TMath::Sign(fSensOffset, x0a);
222     if (x02 * x0a < 0) 
223       x02 = TMath::Sign(fSensOffset, x0a);
224     if (y01 * yhit < 0) 
225       y01 = TMath::Sign(fSensOffset, yhit);
226     if (y02 * yhit < 0) 
227       y02 = TMath::Sign(fSensOffset, yhit);
228     //
229     // find the pads over which the charge distributes
230     GetPadI(x01,y01,fixmin,fiymin);
231     GetPadI(x02,y02,fixmax,fiymax);
232     // 
233     // Set current pad to lower left corner
234     fix=fixmin;
235     fiy=fiymin;
236     SetPadCoord(fix,fiy);
237 }
238
239 void AliMUONSegmentationV1::NextPad()
240 {
241   // 
242   // Step to next pad in integration region
243     if (fix != fixmax) {
244         fix++;
245     } else if (fiy != fiymax) {
246         fix=fixmin;
247         fiy++;
248     } else 
249         printf("\n Error: Stepping outside integration region\n ");
250     SetPadCoord(fix,fiy);
251 }
252
253 Int_t AliMUONSegmentationV1::MorePads()
254 {
255 //
256 // Are there more pads in the integration region
257
258     if (fix == fixmax && fiy == fiymax) {
259         return 0;
260     } else {
261         return 1;       
262     }
263 }
264
265 Int_t AliMUONSegmentationV1::IsParallel2(Int_t iX, Int_t iY)
266 {
267 // test if the pad is read in parallel for zone 2
268 // iX and iY are assumed to be positive and starting at 0 numbering (cF. iX)
269 // returns 1 or 2 if read in parallel, 
270 // according to the actual number in the chain, 0 else
271 //
272 // chainage is        result is 
273 // 1  2  3  1  2  3   1 1 1 2 2 2     y
274 // 7  8  9 10 11 12   0 0 0 0 0 0     ^
275 // 4  5  6  4  5  6   1 1 1 2 2 2     +->x
276 //
277
278     if (iY%3==1) return 0;
279     return (iX%6)/3+1;
280 }
281
282 Int_t AliMUONSegmentationV1::IsParallel3(Int_t iX, Int_t iY)
283 {
284 // test if the pad is read in parallel for zone 3
285 // iX and iY are assumed to be positive and starting at 0 numbering (cF. iX)
286 // returns 1,2 or 3 if read in parallel, 
287 // according to the actual number in the chain, 0 else
288 //
289 // chainage is                 result is
290 //16  2  3  1  2  3  1  2  3   0 1 1 1 2 2 2 3 3
291 // 7  8  9 10 11 12 13 14 15   0 0 0 0 0 0 0 0 0
292 // 4  5  6  4  5  6  4  5  6   1 1 1 2 2 2 3 3 3
293 //
294
295     if (iY%3==1) return 0;
296     return (iX%9)/3+1 - (iY%3==2 && iX%3==0);
297 }
298
299 Int_t AliMUONSegmentationV1::NParallel2(Int_t iX, Int_t iY)
300 {
301 // returns the number of pads connected in parallel for zone 2
302 // iX and iY are assumed to be positive and starting at 0 numbering (cF. iX)
303 //
304 // result is
305 // 2 2 2 2 2 2
306 // 1 1 1 1 1 1
307 // 2 2 2 2 2 2
308 //
309
310     if (iY%3==1) return 1;
311     return 2;
312 }
313
314 Int_t AliMUONSegmentationV1::NParallel3(Int_t iX, Int_t iY)
315 {
316 // test if the pad is read in parallel for zone 3
317 // iX and iY are assumed to be positive and starting at 0 numbering (cF. iX)
318 // returns 1,2 or 3 if read in parallel, 
319 // according to the actual number in the chain, 0 else
320 //
321 // result is 
322 // 1 3 3 2 3 3 2 3 3 
323 // 1 1 1 1 1 1 1 1 1
324 // 3 3 3 3 3 3 3 3 3
325 //
326
327     if (iY%3==1) return 1;
328     if (iY%3==2 && iX%9==0) return 1;
329     return 3 - (iY%3==2 && iX%3==0);
330 }
331
332
333 Int_t AliMUONSegmentationV1::Ix(Int_t trueX, Int_t trueY)
334 {
335 // returns the X number of pad which corresponds to the logical 
336 // channel, expressed in x and y.
337
338     Int_t wix = TMath::Abs(trueX)-1;
339     Int_t wiy = TMath::Abs(trueY)-1;
340     Int_t zone = GetZone(trueX,trueY);
341     Int_t par3;
342     switch (zone) {
343     case 0: return trueX;
344     case 1:
345         if (IsParallel2(wix,wiy) == 2)
346             return (trueX>0)? trueX-3 : trueX+3 ;
347         return trueX;
348     case 2:
349         if ( (par3= IsParallel3(wix,wiy)) )
350             return (trueX>0) ? trueX-3*(par3-1) : trueX+3*(par3-1) ;
351         return trueX ;
352     default :
353         printf("Couille dans AliMUONSegmentationV1::ix\n");
354     }
355     return -1;
356 }
357
358 Int_t AliMUONSegmentationV1::Ix()
359 {
360 // returns the X number of pad which has to increment charge
361 // due to parallel read-out
362 return Ix(fix,fiy);
363 }
364
365 Int_t AliMUONSegmentationV1::ISector()
366 {
367 // This function is of no use for this kind of segmentation.
368     return GetZone(fix,fiy);
369 }
370
371 void AliMUONSegmentationV1::SigGenInit(Float_t x,Float_t y,Float_t z)
372 {
373 //
374 //  Initialises pad and wire position during stepping
375     fxt =x;
376     fyt =y;
377     GetPadI(x,y,fixt,fiyt);
378     fiwt= GetiAnod(x);
379
380 }
381
382 Int_t AliMUONSegmentationV1::SigGenCond(Float_t x,Float_t y,Float_t z)
383 {
384 //
385 //  Signal will be generated if particle crosses pad boundary or
386 //  boundary between two wires. 
387     Int_t ixt;
388     Int_t iyt;
389     GetPadI(x,y,ixt,iyt);
390     Int_t iwt= GetiAnod(x);
391     
392     if ((ixt != fixt) || (iyt !=fiyt) || (iwt != fiwt)) {
393         return 1;
394     } else {
395         return 0;
396     }
397 }
398
399 void AliMUONSegmentationV1::
400 IntegrationLimits(Float_t& x1,Float_t& x2,Float_t& y1, Float_t& y2)
401 {
402 // Get integration limits
403     x1=fxhit-fx-fDpx/2.;
404     x2=x1+fDpx;
405     y1=fyhit-fy-fDpy/2.;
406     y2=y1+fDpy;    
407 }
408
409 void AliMUONSegmentationV1::GetNParallelAndOffset(Int_t iX, Int_t iY,Int_t
410 *Nparallel, Int_t* Offset)
411 {
412 // Get parallel pad
413     Int_t wix = TMath::Abs(iX)-1;
414     Int_t wiy = TMath::Abs(iY)-1;
415     Int_t zone = GetZone(iX,iY);
416     switch (zone) {
417     case 0: 
418         *Nparallel=1;
419         *Offset=0;
420         break;
421     case 1:
422         *Nparallel = NParallel2(wix,wiy);
423         (iX>0) ? *Offset =3 : *Offset = -3;
424         if (IsParallel2(wix,wiy)>1)
425             printf("GetNParallelAndOffset called for existing channel -> answer is crazy\n");
426         break;
427     case 2:
428         *Nparallel = NParallel3(wix,wiy);
429         (iX>0) ? *Offset =3 : *Offset = -3;
430         if (IsParallel3(wix,wiy)>1)
431             printf("GetNParallelAndOffset called for existing channel -> answer is crazy\n");
432         break;
433     }
434 }
435
436
437 Float_t AliMUONSegmentationV1::Distance2AndOffset(Int_t iX, Int_t iY, Float_t X, Float_t Y, Int_t *Offset)
438 {
439 //
440 // Computes the offset for which the physical pad has the minimum distance squared
441 // (returned value) to the given coordinates
442
443     Int_t nPara,offset;
444     GetNParallelAndOffset(iX,iY,&nPara,&offset);
445     Float_t d2min=1E10;
446     for (Int_t i=0;i<nPara; i++)
447     {
448         Float_t x,y;
449         GetPadC(iX+i*offset,iY,x,y);
450         Float_t d2=(x-X)*(x-X) + (y-Y)*(y-Y);
451         if ( d2min > d2)
452         {
453             d2min = d2;
454             *Offset = i*offset;
455         }
456     }
457     return d2min; 
458 }
459
460 void AliMUONSegmentationV1::CleanNeighbours(Int_t* Nlist, Int_t *Xlist, 
461         Int_t *Ylist)
462 {
463 // In the raw neighbours list, some pads do not exist 
464 // and some others are read in parallel ...
465 // So we prune non-existing neighbours from the list (event if this should be
466 // at last not be a problem due to the clustering algorithm...)
467
468     Int_t nTot=0;
469     for (Int_t nList=0;nList<*Nlist;nList++)
470     {
471         // prune if it does not exist
472         if ( Xlist[nList]==0 || Ylist[nList]==0 )
473             continue;
474         // compute true position
475         Xlist[nTot] = Ix(Xlist[nList],Ylist[nList]) ;
476         Ylist[nTot] = Ylist[nList] ;
477         // and prune if it does already exist
478         Int_t nTest;
479         for (nTest=0;nTest<nTot; nTest++)
480         {
481             if ( Xlist[nTest]==Xlist[nTot] && Ylist[nTest]==Ylist[nTot])
482                 // we found it
483                 break ;
484         }
485         if (nTest==nTot)
486             nTot++;
487     }
488     *Nlist = nTot;
489 }
490
491 void AliMUONSegmentationV1::
492 NeighboursNonDiag(Int_t iX, Int_t iY, Int_t* Nlist, Int_t Xlist[12], Int_t Ylist[12])
493 {
494 // returns the X number of pad which has to increment charge
495 // due to parallel read-out
496
497     Int_t nParallel, offset;
498     GetNParallelAndOffset(iX,iY,&nParallel,&offset);
499 //
500 // now fill raw list of neighbours
501     *Nlist=4*nParallel;
502     Xlist[0]=Xlist[1]=iX;Xlist[2]=iX-1;Xlist[3]=iX+1;
503     Ylist[0]=iY-1;Ylist[1]=iY+1;Ylist[2]=Ylist[3]=iY;
504     if (nParallel>1) {
505         Xlist[4]=Xlist[5]=iX+offset;Xlist[6]=iX+offset-1;Xlist[7]=iX+offset+1;
506         Ylist[4]=iY-1;Ylist[5]=iY+1;Ylist[6]=Ylist[7]=iY;
507         if (nParallel>2) {
508             Xlist[8]=Xlist[9]=iX+2*offset;Xlist[10]=iX+2*offset-1;Xlist[11]=iX+2*offset+1;
509             Ylist[8]=iY-1;Ylist[9]=iY+1;Ylist[10]=Ylist[11]=iY;
510         }
511     }
512     CleanNeighbours(Nlist,Xlist,Ylist);
513 }
514
515 void AliMUONSegmentationV1::
516 NeighboursDiag(Int_t iX, Int_t iY, Int_t* Nlist, Int_t Xlist[24], Int_t Ylist[24])
517 {
518 // returns the X number of pad which has to increment charge
519 // due to parallel read-out
520
521     Int_t nParallel, offset;
522     GetNParallelAndOffset(iX,iY,&nParallel,&offset);
523 //
524 // now fill raw list of neighbours
525     *Nlist=0;
526     for (Int_t i=0;i<nParallel;i++)
527         for (Int_t dx=-1;dx<2;dx++)
528             for (Int_t dy=-1;dy<2;dy++)
529             {
530                 if (dx==dy && dy==0)
531                     continue; 
532                 Xlist[*Nlist] = iX + dx + i*offset;
533                 Ylist[*Nlist] = iY + dy;
534                 (*Nlist)++;
535             }
536     CleanNeighbours(Nlist,Xlist,Ylist);
537 }
538
539 void AliMUONSegmentationV1::Neighbours(Int_t iX, Int_t iY, Int_t* Nlist, 
540                                        Int_t Xlist[24], Int_t Ylist[24])
541 {
542 // Get neighbours
543 NeighboursDiag(iX,iY,Nlist,Xlist,Ylist);
544 }
545
546
547 void AliMUONSegmentationV1::GiveTestPoints(Int_t &n, Float_t *x, Float_t *y)
548 {
549 // Return a test point
550     n=1;
551     x[0]=(TMath::Sqrt(frSensMax2)-TMath::Sqrt(frSensMin2))/2/TMath::Sqrt(2.);
552     y[0]=x[0];
553 }
554
555 AliMUONSegmentationV1& AliMUONSegmentationV1::operator =(const AliMUONSegmentationV1 & rhs)
556 {
557 // Dummy assignment operator
558     return *this;
559 }