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