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