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 }