ITS/AliITSgeom.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.4.4.6  2000/06/04 16:33:32  Nilsen
  19 A restructured AliITSgeom class. Now used AliITSgeomMatrix.
  20
  21 Revision 1.4.4.5  2000/03/04 23:42:39  Nilsen
  22 Updated the comments/documentations and improved the maintainability of the
  23 code.
  24
  25 Revision 1.4.4.4  2000/03/02 21:27:07  Nilsen
  26 Added two functions, SetByAngles and SetTrans.
  27
  28 Revision 1.4.4.3  2000/01/23 03:09:10  Nilsen
  29 // fixed compiler warnings for new function LtLErrorMatrix(...)
  30
  31 Revision 1.4.4.2  2000/01/19 23:18:20  Nilsen
  32 Added transformations of Error matrix to AliITSgeom and fixed some typos
  33 in AliITS.h and AliITShitIndex.h
  34
  35 Revision 1.4.4.1  2000/01/12 19:03:32  Nilsen
  36 This is the version of the files after the merging done in December 1999.
  37 See the ReadMe110100.txt file for details
  38
  39 Revision 1.4  1999/10/15 07:03:20  fca
  40 Fixed bug in GetModuleId(Int_t index,Int_t &lay,Int_t &lad, Int_t &det) and
  41 a typo in the creator. aliroot need to be rerun to get a fixed geometry.
  42
  43 Revision 1.3  1999/10/04 15:20:12  fca
  44 Correct syntax accepted by g++ but not standard for static members, remove minor warnings
  45
  46 Revision 1.2  1999/09/29 09:24:20  fca
  47 Introduction of the Copyright and cvs Log
  48
  49 */
  50
  51 ///////////////////////////////////////////////////////////////////////
  52 // ITS geometry manipulation routines.                               //
  53 // Created April 15 1999.                                            //
  54 // version: 0.0.0                                                    //
  55 // By: Bjorn S. Nilsen                                               //
  56 // version: 0.0.1                                                    //
  57 // Updated May 27 1999.                                              //
  58 // Added Cylindrical random and global based changes.               //
  59 // Added  function PrintComparison.                                  //
  60 ///////////////////////////////////////////////////////////////////////
  61
  62
  63 ////////////////////////////////////////////////////////////////////////
  64 //     The local coordinate system by, default, is show in the following
  65 // figures. Also shown are the ladder numbering scheme.
  66 //Begin_Html
  67 /*
  68 <img src="picts/ITS/AliITSgeomMatrix_L1.gif">
  69 </pre>
  70 <br clear=left>
  71 <font size=+2 color=blue>
  72 <p>This shows the relative geometry differences between the ALICE Global
  73 coordinate system and the local detector coordinate system.
  74 </font>
  75 <pre>
  76
  77 <pre>
  78 <img src="picts/ITS/its1+2_convention_front_5.gif">
  79 </pre>
  80 <br clear=left>
  81 <font size=+2 color=blue>
  82 <p>This shows the front view of the SPDs and the orientation of the local
  83 pixel coordinate system. Note that the inner pixel layer has its y coordinate
  84 in the opposite direction from all of the other layers.
  85 </font>
  86 <pre>
  87
  88 <pre>
  89 <img src="picts/ITS/its3+4_convention_front_5.gif">
  90 </pre>
  91 <br clear=left>
  92 <font size=+2 color=blue>
  93 <p>This shows the front view of the SDDs and the orientation of the local
  94 pixel coordinate system.
  95 </font>
  96 <pre>
  97
  98 <pre>
  99 <img src="picts/ITS/its5+6_convention_front_5.gif">
 100 </pre>
 101 <br clear=left>
 102 <font size=+2 color=blue>
 103 <p>This shows the front view of the SSDs and the orientation of the local
 104 pixel coordinate system.
 105 </font>
 106 <pre>
 107 */
 108 //End_Html
 109 //
 110 ////////////////////////////////////////////////////////////////////////
 111
 112 ////////////////////////////////////////////////////////////////////////
 113 //
 114 // version: 0
 115 // Written by Bjorn S. Nilsen
 116 //
 117 // Data Members:
 118 //
 119 // Int_t fNlayers
 120 //     The number of ITS layers for this geometry. By default this
 121 //  is 6, but can be modified by the creator function if there are
 122 // more layers defined.
 123 //
 124 // Int_t *fNlad
 125 //     A pointer to an array fNlayers long containing the number of
 126 // ladders for each layer. This array is typically created and filled
 127 // by the AliITSgeom creator function.
 128 //
 129 // Int_t *fNdet
 130 //     A pointer to an array fNlayers long containing the number of
 131 // active detector volumes for each ladder. This array is typically
 132 // created and filled by the AliITSgeom creator function.
 133 //
 134 // AliITSgeomMatrix *fGm
 135 //     A pointer to an array of AliITSgeomMatrix classes. One element
 136 // per module (detector) in the ITS. AliITSgeomMatrix basicly contains
 137 // all of the necessary information about the detector and it's coordinate
 138 // transformations.
 139 //
 140 // TObjArray *fShape
 141 //     A pointer to an array of TObjects containing the detailed shape
 142 // information for each type of detector used in the ITS. For example
 143 // I have created AliITSgeomSPD, AliITSgeomSDD, and AliITSgeomSSD as
 144 // example structures, derived from TObjects, to hold the detector
 145 // information. I would recommend that one element in each of these
 146 // structures, that which describes the shape of the active volume,
 147 // be one of the ROOT classes derived from TShape. In this way it would
 148 // be easy to have the display program display the correct active
 149 // ITS volumes. See the example classes AliITSgeomSPD, AliITSgeomSDD,
 150 // and AliITSgeomSSD for a more detailed example.
 151 ////////////////////////////////////////////////////////////////////////
 152 #include <iostream.h>
 153 #include <fstream.h>
 154 #include <iomanip.h>
 155 #include <stdio.h>
 156 #include <string.h>
 157 #include <ctype.h>
 158 #include <TObject.h>
 159 #include <TRandom.h>
 160
 161 #include "AliITSgeom.h"
 162 #include "AliITSgeomMatrix.h"
 163 #include "AliITSgeomSPD.h"
 164 #include "AliITSgeomSDD.h"
 165 #include "AliITSgeomSSD.h"
 166
 167 ClassImp(AliITSgeom)
 168
 169 //_____________________________________________________________________
 170 AliITSgeom::AliITSgeom(){
 171 ////////////////////////////////////////////////////////////////////////
 172 //     The default constructor for the AliITSgeom class. It, by default,
 173 // sets fNlayers to zero and zeros all pointers.
 174 ////////////////////////////////////////////////////////////////////////
 175   // Default constructor.
 176   // Do not allocate anything zero everything
 177    fTrans   = 0; // standard GEANT global/local coordinate system.
 178    fNlayers = 0;
 179    fNlad    = 0;
 180    fNdet    = 0;
 181    fGm      = 0;
 182    fShape   = 0;
 183    return;
 184 }
 185
 186 //_____________________________________________________________________
 187 AliITSgeom::~AliITSgeom(){
 188 ////////////////////////////////////////////////////////////////////////
 189 //     The destructor for the AliITSgeom class. If the arrays fNlad,
 190 // fNdet, or fGm have had memory allocated to them, there pointer values
 191 // are non zero, then this memory space is freed and they are set
 192 // to zero. In addition, fNlayers is set to zero. The destruction of
 193 // TObjArray fShape is, by default, handled by the TObjArray destructor.
 194 ////////////////////////////////////////////////////////////////////////
 195   // Default destructor.
 196   // if arrays exist delete them. Then set everything to zero.
 197    if(fGm!=0){
 198       for(Int_t i=0;i<fNlayers;i++) delete fGm[i];
 199       delete fGm;
 200    } // end if fGm!=0
 201    if(fNlad!=0) delete[] fNlad;
 202    if(fNdet!=0) delete[] fNdet;
 203    fNlayers = 0;
 204    fNlad    = 0;
 205    fNdet    = 0;
 206    fGm      = 0;
 207    return;
 208 }
 209 //______________________________________________________________________
 210 void AliITSgeom::ReadNewFile(const char *filename){
 211     printf("New file format not defined yet\n");
 212     return;
 213 }
 214 //_____________________________________________________________________
 215 AliITSgeom::AliITSgeom(const char *filename){
 216 ////////////////////////////////////////////////////////////////////////
 217 //     The constructor for the AliITSgeom class. All of the data to fill
 218 // this structure is read in from the file given my the input filename.
 219 ////////////////////////////////////////////////////////////////////////
 220    FILE     *pf;
 221    Int_t    i,lm=0,id[3];
 222    Int_t    l,a,d;
 223    Float_t  x,y,z,o,p,q,r,s,t;
 224    Double_t rot6[6],tran[3];
 225    char     buf[200],*buff=0; // input character buffer;
 226
 227    pf = fopen(filename,"r");
 228
 229    fNlayers = 6; // set default number of ladders
 230 TryAgain:
 231    fNlad    = new Int_t[fNlayers];
 232    fNdet    = new Int_t[fNlayers];
 233    fNmodules = 0;
 234    // find the number of ladders and detectors in this geometry.
 235    for(i=0;i<fNlayers;i++){fNlad[i]=fNdet[i]=0;} // zero out arrays
 236    while(fgets(buf,200,pf)!=NULL){ // for ever loop
 237       for(i=0;i<200;i++)if(buf[i]!=' '){ // remove blank spaces.
 238            buff = &(buf[i]);
 239            break;
 240       } // end for i
 241       // remove blank lines and comments.
 242       if(buff[0]=='\n'||buff[0]=='#'||buff[0]=='!'||
 243          (buff[0]=='/'&&buff[1]=='/')) continue;
 244       if(isalpha(buff[0])) { // must be the new file formated file.
 245             fclose(pf);
 246             delete[] fNlad;delete[] fNdet;
 247             ReadNewFile(filename);
 248             return;
 249       } // end if isalpha(buff[0])
 250       sscanf(buff,"%d %d %d %f %f %f %f %f %f %f %f %f",
 251                   &l,&a,&d,&x,&y,&z,&o,&p,&q,&r,&s,&t);
 252       if(l>lm) lm = l;
 253       if(l<1 || l>fNlayers) {
 254          printf("error in file %s layer=%d min. is 1 max is %d/n",
 255                  filename,l,fNlayers);
 256          continue;
 257       }// end if l
 258       fNmodules++;
 259       if(l<=fNlayers&&fNlad[l-1]<a) fNlad[l-1] = a;
 260       if(l<=fNlayers&&fNdet[l-1]<d) fNdet[l-1] = d;
 261    } // end while ever loop
 262    if(lm>fNlayers){
 263         delete[] fNlad;
 264         delete[] fNdet;
 265         fNlayers = lm;
 266         goto TryAgain;
 267    } // end if lm>fNlayers
 268    // counted the number of ladders and detectors now allocate space.
 269    fGm = new AliITSgeomMatrix*[fNmodules];
 270
 271    // Set up Shapes for a default configuration of 6 layers.
 272    fTrans   = 0; // standard GEANT global/local coordinate system.
 273    fShape = new TObjArray(3);
 274    AddShape((TObject *) new AliITSgeomSPD());  // shape 0
 275    AddShape((TObject *) new AliITSgeomSDD());  // shape 1
 276    AddShape((TObject *) new AliITSgeomSPD());  // shape 2
 277
 278    // prepare to read in transforms
 279    lm = 0; // reuse lm as counter of modules.
 280    rewind(pf); // start over reading file
 281    while(fgets(buf,200,pf)!=NULL){ // for ever loop
 282       for(i=0;i<200;i++)if(buf[i]!=' '){ // remove blank spaces.
 283            buff = &(buf[i]);
 284            break;
 285       } // end for i
 286       // remove blank lines and comments.
 287       if(buff[0]=='\n'||buff[0]=='#'||buff[0]=='!'||
 288         (buff[0]=='/'&&buff[1]=='/')) continue;
 289       x = y = z = o = p = q = r = s = t = 0.0;
 290       sscanf(buff,"%d %d %d %f %f %f %f %f %f %f %f %f",
 291                   &l,&a,&d,&x,&y,&z,&o,&p,&q,&r,&s,&t);
 292       if(l<1 || l>fNlayers) {
 293          printf("error in file %s layer=%d min. is 1 max is %d/n",
 294                  filename,l,fNlayers);
 295          continue;
 296       }// end if l
 297       id[0] = l;id[1] = a;id[2] = d;
 298       tran[0] = tran[1] = tran[2]  = 0.0;
 299       tran[0] = (Double_t)x;tran[1] = (Double_t)y;tran[2] = (Double_t)z;
 300       rot6[0] = rot6[1] = rot6[2] = rot6[3] = rot6[4] = rot6[5] =0.0;
 301       rot6[0] = (Double_t)o;rot6[1] = (Double_t)p;rot6[2] = (Double_t)q;
 302       rot6[3] = (Double_t)r;rot6[4] = (Double_t)s;rot6[5] = (Double_t)t;
 303       switch (l){
 304       case 1: case 2: // layer 1 or2 SPD
 305           fGm[lm++] = new AliITSgeomMatrix(rot6,0,id,tran);
 306           break;
 307       case 3: case 4: // layer 3 or 4 SDD
 308           fGm[lm++] = new AliITSgeomMatrix(rot6,1,id,tran);
 309           break;
 310       case 5: case 6: // layer 5 or 6 SSD
 311           fGm[lm++] = new AliITSgeomMatrix(rot6,2,id,tran);
 312           break;
 313       } // end switch
 314    } // end while ever loop
 315    fclose(pf);
 316 }
 317
 318 //________________________________________________________________________
 319 AliITSgeom::AliITSgeom(AliITSgeom &source){
 320 ////////////////////////////////////////////////////////////////////////
 321 //     The copy constructor for the AliITSgeom class. It calls the
 322 // = operator function. See the = operator function for more details.
 323 ////////////////////////////////////////////////////////////////////////
 324
 325     *this = source;  // Just use the = operator for now.
 326
 327     return;
 328 }
 329
 330 //________________________________________________________________________
 331 void AliITSgeom::operator=(AliITSgeom &source){
 332 ////////////////////////////////////////////////////////////////////////
 333 //     The = operator function for the AliITSgeom class. It makes an
 334 // independent copy of the class in such a way that any changes made
 335 // to the copied class will not affect the source class in any way.
 336 // This is required for many ITS alignment studies where the copied
 337 // class is then modified by introducing some misalignment.
 338 ////////////////////////////////////////////////////////////////////////
 339    Int_t i;
 340
 341    if(this == &source) return; // don't assign to ones self.
 342
 343    // if there is an old structure allocated delete it first.
 344    if(this->fGm != 0){
 345       for(i=0;i<this->fNmodules;i++) delete this->fGm[i];
 346       delete this->fGm;
 347    } // end if fGm != 0
 348    if(fNlad != 0) delete[] fNlad;
 349    if(fNdet != 0) delete[] fNdet;
 350
 351    this->fTrans    = source.fTrans;
 352    this->fNmodules = source.fNmodules;
 353    this->fNlayers = source.fNlayers;
 354    this->fNlad = new Int_t[fNlayers];
 355    for(i=0;i<this->fNlayers;i++) this->fNlad[i] = source.fNlad[i];
 356    this->fNdet = new Int_t[fNlayers];
 357    for(i=0;i<this->fNlayers;i++) this->fNdet[i] = source.fNdet[i];
 358    this->fShape = new TObjArray(*(source.fShape));//This does not make a proper copy.
 359    this->fGm = new AliITSgeomMatrix*[this->fNmodules];
 360    for(i=0;i<this->fNmodules;i++){
 361        this->fGm[i] = new AliITSgeomMatrix(*(source.fGm[i]));
 362    } // end for i
 363    return;
 364 }//_____________________________________________________________________
 365 Int_t AliITSgeom::GetModuleIndex(const Int_t lay,const Int_t lad,
 366                                  const Int_t det){
 367 ////////////////////////////////////////////////////////////////////////
 368 //      This routine computes the module index number from the layer,
 369 // ladder, and detector numbers. The number of ladders and detectors
 370 // per layer is determined when this geometry package is constructed,
 371 // see AliITSgeom(const char *filename) for specifics.
 372 ////////////////////////////////////////////////////////////////////////
 373     Int_t i,j,k,id[3];
 374
 375     i = fNdet[lay-1] * (lad-1) + det - 1;
 376     j = 0;
 377     for(k=0;k<lay-1;k++) j += fNdet[k]*fNlad[k];
 378     i = i+j;
 379     fGm[i]->GetIndex(id);
 380     if(id[0]==lay&&id[1]==lad&&id[2]==det) return i;
 381     // Array of modules fGm is not in expected order. Search for this index
 382     for(i=0;i<fNmodules;i++){
 383         fGm[i]->GetIndex(id);
 384         if(id[0]==lay&&id[1]==lad&&id[2]==det) return i;
 385     } // end for i
 386     // This layer ladder and detector combination does not exist return -1.
 387     return -1;
 388 }
 389 //______________________________________________________________________
 390 void AliITSgeom::GetModuleId(const Int_t index,
 391                              Int_t &lay,Int_t &lad,Int_t &det){
 392 ////////////////////////////////////////////////////////////////////////
 393 //      This routine computes the layer, ladder and detector number
 394 // given the module index number. The number of ladders and detectors
 395 // per layer is determined when this geometry package is constructed,
 396 // see AliITSgeom(const char *filename) for specifics.
 397 ////////////////////////////////////////////////////////////////////////
 398     Int_t id[3];
 399
 400     fGm[index]->GetIndex(id);
 401     lay = id[0]; lad = id[1]; det = id[2];
 402     return;
 403
 404     // The old way kept for posterity.
 405 /*
 406     Int_t i,j,k;
 407     j = 0;
 408     for(k=0;k<fNlayers;k++){
 409         j += fNdet[k]*fNlad[k];
 410         if(j>index)break;
 411     } // end for k
 412     lay = k+1;
 413     i = index -j + fNdet[k]*fNlad[k];
 414     j = 0;
 415     for(k=0;k<fNlad[lay-1];k++){
 416         j += fNdet[lay-1];
 417         if(j>i)break;
 418     } // end for k
 419     lad = k+1;
 420     det = 1+i-fNdet[lay-1]*k;
 421     return;
 422 */
 423 }
 424 //___________________________________________________________________________
 425 Int_t AliITSgeom::GetStartDet(const Int_t dtype){
 426   /////////////////////////////////////////////////////////////////////////
 427   // returns the starting module index value for a give type of detector id
 428   /////////////////////////////////////////////////////////////////////////
 429
 430   switch(dtype){
 431   case 0:
 432      return GetModuleIndex(1,1,1);
 433      break;
 434   case 1:
 435      return GetModuleIndex(3,1,1);
 436      break;
 437   case 2:
 438      return GetModuleIndex(5,1,1);
 439      break;
 440   default:
 441      printf("<AliITSgeom::GetFirstDet> undefined detector type\n");
 442      return 0;
 443   } // end switch
 444
 445   printf("<AliITSgeom::GetFirstDet> undefined detector type\n");
 446   return 0;
 447 }
 448
 449 //___________________________________________________________________________
 450 Int_t AliITSgeom::GetLastDet(const Int_t dtype){
 451   /////////////////////////////////////////////////////////////////////////
 452   // returns the last module index value for a give type of detector id
 453   /////////////////////////////////////////////////////////////////////////
 454
 455   switch(dtype){
 456   case 0:
 457      return GetLastSPD();
 458      break;
 459    case 1:
 460      return GetLastSDD();
 461      break;
 462    case 2:
 463      return GetLastSSD();
 464      break;
 465    default:
 466      printf("<AliITSgeom::GetLastDet> undefined detector type\n");
 467      return 0;
 468   } // end switch
 469
 470   printf("<AliITSgeom::GetLastDet> undefined detector type\n");
 471   return 0;
 472 }
 473
 474 //___________________________________________________________________________
 475 void AliITSgeom::PrintComparison(FILE *fp,AliITSgeom *other){
 476 ////////////////////////////////////////////////////////////////////////
 477 //     This function was primarily created for diagnostic reasons. It
 478 // print to a file pointed to by the file pointer fp the difference
 479 // between two AliITSgeom classes. The format of the file is basicly,
 480 // define d? to be the difference between the same element of the two
 481 // classes. For example dfrx = this->fGm[i][j].frx - other->fGm[i][j].frx.
 482 // if(at least one of dfx0, dfy0, dfz0,dfrx,dfry,dfrz are non zero) then print
 483 // layer ladder detector dfx0 dfy0 dfz0 dfrx dfry dfrz
 484 // if(at least one of the 9 elements of dfr[] are non zero) then print
 485 // layer ladder detector dfr[0] dfr[1] dfr[2]
 486 //                       dfr[3] dfr[4] dfr[5]
 487 //                       dfr[6] dfr[7] dfr[8]
 488 // Only non zero values are printed to save space. The differences are
 489 // typical written to a file because there are usually a lot of numbers
 490 // printed out and it is usually easier to read them in some nice editor
 491 // rather than zooming quickly past you on a screen. fprintf is used to
 492 // do the printing. The fShapeIndex difference is not printed at this time.
 493 ////////////////////////////////////////////////////////////////////////
 494    Int_t    i,j,idt[3],ido[3];
 495    Double_t tt[3],to[3];  // translation
 496    Double_t rt[3],ro[3];  // phi in radians
 497    Double_t mt[3][3],mo[3][3]; // matrixes
 498    AliITSgeomMatrix *gt,*go;
 499    Bool_t   t;
 500
 501    for(i=0;i<this->fNmodules;i++){
 502          gt  =  this->GetGeomMatrix(i);
 503          go  = other->GetGeomMatrix(i);
 504          gt->GetIndex(idt);
 505          go->GetIndex(ido);
 506          t = kFALSE;
 507          for(i=0;i<3;i++) t = t&&idt[i]!=ido[i];
 508          if(t) fprintf(fp,"%4.4d %1.1d %2.2d %2.2d %1.1d %2.2d %2.2d\n",i,
 509                        idt[0],idt[1],idt[2],ido[0],ido[1],ido[2]);
 510          gt->GetTranslation(tt);
 511          go->GetTranslation(to);
 512          gt->GetAngles(rt);
 513          go->GetAngles(ro);
 514          t = kFALSE;
 515          for(i=0;i<3;i++) t = t&&tt[i]!=to[i];
 516          if(t) fprintf(fp,"%1.1d %2.2d %2.2d dTrans=%f %f %f drot=%f %f %f\n",
 517                        idt[0],idt[1],idt[2],
 518                        tt[0]-to[0],tt[1]-to[1],tt[2]-to[2],
 519                        rt[0]-ro[0],rt[1]-ro[1],rt[2]-ro[2]);
 520          t = kFALSE;
 521          gt->GetMatrix(mt);
 522          go->GetMatrix(mo);
 523          for(i=0;i<3;i++)for(j=0;j<3;j++)  t = mt[i][j] != mo[i][j];
 524          if(t){
 525              fprintf(fp,"%1.1d %2.2d %2.2d dfr= %e %e %e\n",
 526                      idt[0],idt[1],idt[2],
 527                  mt[0][0]-mo[0][0],mt[0][1]-mo[0][1],mt[0][2]-mo[0][2]);
 528              fprintf(fp,"        dfr= %e %e %e\n",
 529                  mt[1][0]-mo[1][0],mt[1][1]-mo[1][1],mt[1][2]-mo[1][2]);
 530              fprintf(fp,"        dfr= %e %e %e\n",
 531                  mt[2][0]-mo[2][0],mt[2][1]-mo[2][1],mt[2][2]-mo[2][2]);
 532          } // end if t
 533    } // end for i
 534    return;
 535 }
 536
 537 //___________________________________________________________________________
 538 void AliITSgeom::PrintData(FILE *fp,
 539                            const Int_t lay,const Int_t lad,const Int_t det){
 540 ////////////////////////////////////////////////////////////////////////
 541 //     This function prints out the coordinate transformations for
 542 // the particular detector defined by layer, ladder, and detector
 543 // to the file pointed to by the File pointer fp. fprintf statements
 544 // are used to print out the numbers. The format is
 545 // layer ladder detector Trans= fx0 fy0 fz0 rot= frx fry frz Shape=fShapeIndex
 546 //                         dfr= fr[0] fr[1] fr[2]
 547 //                         dfr= fr[3] fr[4] fr[5]
 548 //                         dfr= fr[6] fr[7] fr[8]
 549 // By indicating which detector, some control over the information
 550 // is given to the user. The output it written to the file pointed
 551 // to by the file pointer fp. This can be set to stdout if you want.
 552 ////////////////////////////////////////////////////////////////////////
 553    AliITSgeomMatrix *gt;
 554    Double_t t[3],r[3],m[3][3];
 555
 556    gt = this->GetGeomMatrix(GetModuleIndex(lay,lad,det));
 557    gt->GetTranslation(t);
 558    gt->GetAngles(r);
 559    fprintf(fp,"%1.1d %2.2d %2.2d Trans=%f %f %f rot=%f %f %f Shape=%d\n",
 560            lay,lad,det,t[0],t[1],t[2],r[0],r[1],r[2],
 561            gt->GetDetectorIndex());
 562    gt->GetMatrix(m);
 563    fprintf(fp,"        dfr= %e %e %e\n",m[0][0],m[0][1],m[0][2]);
 564    fprintf(fp,"        dfr= %e %e %e\n",m[1][0],m[1][1],m[1][2]);
 565    fprintf(fp,"        dfr= %e %e %e\n",m[2][0],m[2][1],m[2][2]);
 566    return;
 567 }
 568 //___________________________________________________________________________
 569 ofstream & AliITSgeom::PrintGeom(ofstream &R__b){
 570 ////////////////////////////////////////////////////////////////////////
 571 //     The default Streamer function "written by ROOT" doesn't write out
 572 // the arrays referenced by pointers. Therefore, a specific Streamer function
 573 // has to be written. This function should not be modified but instead added
 574 // on to so that older versions can still be read. The proper handling of
 575 // the version dependent streamer function hasn't been written do to the lack
 576 // of finding an example at the time of writing.
 577 ////////////////////////////////////////////////////////////////////////
 578    // Stream an object of class AliITSgeom.
 579     Int_t i;
 580
 581     R__b.setf(ios::scientific);
 582     R__b << fTrans << " ";
 583     R__b << fNmodules << " ";
 584     R__b << fNlayers << " ";
 585     for(i=0;i<fNlayers;i++) R__b << fNlad[i] << " ";
 586     for(i=0;i<fNlayers;i++) R__b << fNdet[i] << "\n";
 587     for(i=0;i<fNmodules;i++) {
 588         R__b <<setprecision(16) << *(fGm[i]) << "\n";
 589     } // end for i
 590     return R__b;
 591 }
 592 //___________________________________________________________________________
 593 ifstream & AliITSgeom::ReadGeom(ifstream &R__b){
 594 ////////////////////////////////////////////////////////////////////////
 595 //     The default Streamer function "written by ROOT" doesn't write out
 596 // the arrays referenced by pointers. Therefore, a specific Streamer function
 597 // has to be written. This function should not be modified but instead added
 598 // on to so that older versions can still be read. The proper handling of
 599 // the version dependent streamer function hasn't been written do to the lack
 600 // of finding an example at the time of writing.
 601 ////////////////////////////////////////////////////////////////////////
 602    // Stream an object of class AliITSgeom.
 603       Int_t i;
 604
 605       fNlad = new Int_t[fNlayers];
 606       fNdet = new Int_t[fNlayers];
 607       if(fGm!=0){
 608           for(i=0;i<fNmodules;i++) delete fGm[i];
 609           delete fGm;
 610       } // end if fGm!=0
 611
 612       R__b >> fTrans >> fNmodules >> fNlayers;
 613       fNlad = new Int_t[fNlayers];
 614       fNdet = new Int_t[fNlayers];
 615       for(i=0;i<fNlayers;i++) R__b >> fNlad[i];
 616       for(i=0;i<fNlayers;i++) R__b >> fNdet[i];
 617       fGm = new AliITSgeomMatrix*[fNmodules];
 618       for(i=0;i<fNmodules;i++){
 619           fGm[i] = new AliITSgeomMatrix;
 620           R__b >> *(fGm[i]);
 621       } // end for i
 622       return R__b;
 623 }
 624 //___________________________________________________________________________
 625 void AliITSgeom::Streamer(TBuffer &R__b){
 626 ////////////////////////////////////////////////////////////////////////
 627 //     The default Streamer function "written by ROOT" doesn't write out
 628 // the arrays referenced by pointers. Therefore, a specific Streamer function
 629 // has to be written. This function should not be modified but instead added
 630 // on to so that older versions can still be read. The proper handling of
 631 // the version dependent streamer function hasn't been written do to the lack
 632 // of finding an example at the time of writing.
 633 ////////////////////////////////////////////////////////////////////////
 634     // Stream an object of class AliITSgeom.
 635     Int_t  i,j,k,l;
 636     UInt_t R__s=0, R__c=0;
 637
 638     if (R__b.IsReading()) {
 639         Version_t R__v = R__b.ReadVersion();
 640         if (R__v==1) {
 641             if(fNlad!=0) delete[] fNlad;
 642             if(fNdet!=0) delete[] fNdet;
 643             if(fGm!=0){
 644                 for(i=0;i<fNlayers;i++) delete[] fGm[i];
 645                 delete[] fGm;
 646             } // end if fGm!=0
 647             Int_t    idt,id[3],inmax;
 648             Double_t t[3],r[3],m[9],s[3][3];
 649             TObject::Streamer(R__b);
 650             fTrans = 0;
 651             R__b >> fNlayers;
 652             fNlad = new Int_t[fNlayers];
 653             fNdet = new Int_t[fNlayers];
 654             for(i=0;i<fNlayers;i++) R__b >> fNlad[i];
 655             for(i=0;i<fNlayers;i++) R__b >> fNdet[i];
 656             fNmodules = GetModuleIndex(fNlayers,fNlad[fNlayers-1],
 657                                        fNdet[fNlayers-1]);
 658             fGm = new AliITSgeomMatrix*[fNmodules];
 659             inmax = 0;
 660             for(i=0;i<fNlayers;i++){
 661                 for(j=0;j<fNlad[i]*fNdet[i];j++){
 662                     R__b >> idt;
 663                     R__b >> t[0];
 664                     R__b >> t[1];
 665                     R__b >> t[2];
 666                     R__b >> r[0];
 667                     R__b >> r[1];
 668                     R__b >> r[2];
 669                     for(k=0;k<9;k++) R__b >> m[k];
 670                     for(k=0;k<3;k++)for(l=0;l<3;l++) s[k][l] = m[3*k+l];
 671                     GetModuleId(inmax,id[0],id[1],id[2]);
 672                     fGm[inmax++] = new AliITSgeomMatrix(idt,id,s,t);
 673                 } // end for j
 674             } // end for i
 675             R__b >> fShape;
 676         } else if(R__v==2){
 677             if(fNlad!=0) delete[] fNlad;
 678             if(fNdet!=0) delete[] fNdet;
 679             if(fGm!=0){for(i=0;i<fNmodules;i++) delete fGm[i];delete[] fGm;}
 680             TObject::Streamer(R__b);
 681             R__b >> fTrans;
 682             R__b >> fNlayers;
 683             R__b.ReadArray(fNlad);
 684             R__b.ReadArray(fNdet);
 685             R__b >> fShape;
 686             R__b >> fNmodules;
 687             fGm = new AliITSgeomMatrix*[fNmodules];
 688             for(i=0;i<fNmodules;i++){
 689                 fGm[i] = new AliITSgeomMatrix;
 690                 fGm[i]->Streamer(R__b);
 691             } // end for i
 692             //R__b.ReadArray(fGm);
 693             R__b.CheckByteCount(R__s, R__c, AliITSgeom::IsA());
 694         } // end if R__v==?
 695     } else { // Writing.
 696         R__c = R__b.WriteVersion(AliITSgeom::IsA(), kTRUE);
 697         TObject::Streamer(R__b);
 698         R__b << fTrans;
 699         R__b << fNlayers;
 700         R__b.WriteArray(fNlad, fNlayers);
 701         R__b.WriteArray(fNdet, fNlayers);
 702         R__b << fShape;
 703         R__b << fNmodules;
 704         //R__b.WriteArray(fGm, __COUNTER__);
 705         for(i=0;i<fNmodules;i++){
 706             fGm[i]->Streamer(R__b);
 707         } // end for i
 708         R__b.SetByteCount(R__c, kTRUE);
 709     } // end if reading/writing.
 710 }
 711 //______________________________________________________________________
 712 //     The following routines modify the transformation of "this"
 713 // geometry transformations in a number of different ways.
 714 //______________________________________________________________________
 715 void AliITSgeom::GlobalChange(const Float_t *tran,const Float_t *rot){
 716 ////////////////////////////////////////////////////////////////////////
 717 //     This function performs a Cartesian translation and rotation of
 718 // the full ITS from its default position by an amount determined by
 719 // the three element arrays dtranslation and drotation. If every element
 720 // of dtranslation and drotation are zero then there is no change made
 721 // the geometry. The change is global in that the exact same translation
 722 // and rotation is done to every detector element in the exact same way.
 723 // The units of the translation are those of the Monte Carlo, usually cm,
 724 // and those of the rotation are in radians. The elements of dtranslation
 725 // are dtranslation[0] = x, dtranslation[1] = y, and dtranslation[2] = z.
 726 // The elements of drotation are drotation[0] = rx, drotation[1] = ry, and
 727 // drotation[2] = rz. A change in x will move the hole ITS in the ALICE
 728 // global x direction, the same for a change in y. A change in z will
 729 // result in a translation of the ITS as a hole up or down the beam line.
 730 // A change in the angles will result in the inclination of the ITS with
 731 // respect to the beam line, except for an effective rotation about the
 732 // beam axis which will just rotate the ITS as a hole about the beam axis.
 733 ////////////////////////////////////////////////////////////////////////
 734    Int_t    i,j;
 735    Double_t t[3],r[3];
 736    AliITSgeomMatrix *g;
 737
 738    fTrans = (fTrans && 0xfffd) + 2;  // set bit 1 true.
 739    for(i=0;i<fNmodules;i++){
 740          g = this->GetGeomMatrix(i);
 741          g->GetTranslation(t);
 742          g->GetAngles(r);
 743          for(j=0;j<3;j++){
 744               t[j] += tran[j];
 745               r[j] += rot[j];
 746          } // end for j
 747          g->SetTranslation(t);
 748          g->SetAngles(r);
 749    } // end for i
 750    return;
 751 }
 752 //___________________________________________________________________________
 753 void AliITSgeom::GlobalCylindericalChange(const Float_t *tran,const Float_t *rot){
 754 ////////////////////////////////////////////////////////////////////////
 755 //     This function performs a cylindrical translation and rotation of
 756 // each ITS element by a fixed about in radius, rphi, and z from its
 757 // default position by an amount determined by the three element arrays
 758 // dtranslation and drotation. If every element of dtranslation and
 759 // drotation are zero then there is no change made the geometry. The
 760 // change is global in that the exact same distance change in translation
 761 // and rotation is done to every detector element in the exact same way.
 762 // The units of the translation are those of the Monte Carlo, usually cm,
 763 // and those of the rotation are in radians. The elements of dtranslation
 764 // are dtranslation[0] = r, dtranslation[1] = rphi, and dtranslation[2] = z.
 765 // The elements of drotation are drotation[0] = rx, drotation[1] = ry, and
 766 // drotation[2] = rz. A change in r will results in the increase of the
 767 // radius of each layer by the same about. A change in rphi will results in
 768 // the rotation of each layer by a different angle but by the same
 769 // circumferential distance. A change in z will result in a translation
 770 // of the ITS as a hole up or down the beam line. A change in the angles
 771 // will result in the inclination of the ITS with respect to the beam
 772 // line, except for an effective rotation about the beam axis which will
 773 // just rotate the ITS as a hole about the beam axis.
 774 ////////////////////////////////////////////////////////////////////////
 775    Int_t    i,j;
 776    Double_t t[3],ro[3],r,r0,phi,rphi;
 777    AliITSgeomMatrix *g;
 778
 779    fTrans = (fTrans && 0xfffd) + 2;  // set bit 1 true.
 780    for(i=0;i<fNmodules;i++){
 781          g = this->GetGeomMatrix(i);
 782          g->GetTranslation(t);
 783          g->GetAngles(ro);
 784          r = r0= TMath::Hypot(t[1],t[0]);
 785          phi   = TMath::ATan2(t[1],t[0]);
 786          rphi  = r0*phi;
 787          r    += tran[0];
 788          rphi += tran[1];
 789          phi   = rphi/r0;
 790          t[0]  = r*TMath::Cos(phi);
 791          t[1]  = r*TMath::Sin(phi);
 792          t[2] += tran[2];
 793          for(j=0;j<3;j++){
 794               ro[j] += rot[j];
 795          } // end for j
 796          g->SetTranslation(t);
 797          g->SetAngles(ro);
 798    } // end for i
 799    return;
 800 }
 801 //___________________________________________________________________________
 802 void AliITSgeom::RandomChange(const Float_t *stran,const Float_t *srot){
 803 ////////////////////////////////////////////////////////////////////////
 804 //     This function performs a Gaussian random displacement and/or
 805 // rotation about the present global position of each active
 806 // volume/detector of the ITS. The sigma of the random displacement
 807 // is determined by the three element array stran, for the
 808 // x y and z translations, and the three element array srot,
 809 // for the three rotation about the axis x y and z.
 810 ////////////////////////////////////////////////////////////////////////
 811    Int_t    i,j;
 812    Double_t t[3],r[3];
 813    TRandom ran;
 814    AliITSgeomMatrix *g;
 815
 816    fTrans = (fTrans && 0xfffd) + 2;  // set bit 1 true.
 817    for(i=0;i<fNmodules;i++){
 818          g = this->GetGeomMatrix(i);
 819          g->GetTranslation(t);
 820          g->GetAngles(r);
 821          for(j=0;j<3;j++){
 822               t[j] += ran.Gaus(0.0,stran[j]);
 823               r[j] += ran.Gaus(0.0, srot[j]);
 824          } // end for j
 825          g->SetTranslation(t);
 826          g->SetAngles(r);
 827    } // end for i
 828    return;
 829 }
 830 //___________________________________________________________________________
 831 void AliITSgeom::RandomCylindericalChange(const Float_t *stran,
 832                                           const Float_t *srot){
 833 ////////////////////////////////////////////////////////////////////////
 834 //     This function performs a Gaussian random displacement and/or
 835 // rotation about the present global position of each active
 836 // volume/detector of the ITS. The sigma of the random displacement
 837 // is determined by the three element array stran, for the
 838 // r rphi and z translations, and the three element array srot,
 839 // for the three rotation about the axis x y and z. This random change
 840 // in detector position allow for the simulation of a random uncertainty
 841 // in the detector positions of the ITS.
 842 ////////////////////////////////////////////////////////////////////////
 843    Int_t    i,j;
 844    Double_t t[3],ro[3],r,r0,phi,rphi;
 845    TRandom ran;
 846    AliITSgeomMatrix *g;
 847
 848    fTrans = (fTrans && 0xfffd) + 2;  // set bit 1 true.
 849    for(i=0;i<fNmodules;i++){
 850          g = this->GetGeomMatrix(i);
 851          g->GetTranslation(t);
 852          g->GetAngles(ro);
 853          r = r0= TMath::Hypot(t[1],t[0]);
 854          phi   = TMath::ATan2(t[1],t[0]);
 855          rphi  = r0*phi;
 856          r    += ran.Gaus(0.0,stran[0]);
 857          rphi += ran.Gaus(0.0,stran[1]);
 858          phi   = rphi/r0;
 859          t[0]  = r*TMath::Cos(phi);
 860          t[1]  = r*TMath::Sin(phi);
 861          t[2] += ran.Gaus(0.0,stran[2]);
 862          for(j=0;j<3;j++){
 863               ro[j] += ran.Gaus(0.0, srot[j]);
 864          } // end for j
 865          g->SetTranslation(t);
 866          g->SetAngles(ro);
 867    } // end for i
 868    return;
 869 }
 870 //______________________________________________________________________
 871 void AliITSgeom::GeantToTracking(AliITSgeom &source){
 872 /////////////////////////////////////////////////////////////////////////
 873 //     Copy the geometry data but change it to go between the ALICE
 874 // Global coordinate system to that used by the ITS tracking. A slightly
 875 // different coordinate system is used when tracking. This coordinate
 876 // system is only relevant when the geometry represents the cylindrical
 877 // ALICE ITS geometry. For tracking the Z axis is left alone but X-> -Y
 878 // and Y-> X such that X always points out of the ITS cylinder for every
 879 // layer including layer 1 (where the detectors are mounted upside down).
 880 //Begin_Html
 881 /*
 882 <img src="picts/ITS/AliITSgeomMatrix_T1.gif">
 883 */
 884 //End_Html
 885 ////////////////////////////////////////////////////////////////////////
 886    Int_t    i,j,k,l,id[3];
 887    Double_t R0[3][3],R1[3][3];
 888    Double_t A0[3][3] = {{0.,+1.,0.},{-1.,0.,0.},{0.,0.,+1.}};
 889    Double_t A1[3][3] = {{0.,-1.,0.},{+1.,0.,0.},{0.,0.,+1.}};
 890
 891    *this = source;  // copy everything
 892    for(i=0;i<GetIndexMax();i++){
 893        fGm[i]->GetIndex(id);
 894        fGm[i]->GetMatrix(R0);
 895        if(id[0]==1){ // Layer 1 is treated different from the others.
 896            for(j=0;j<3;j++) for(k=0;k<3;k++){
 897                R1[j][k] = 0.;
 898                for(l=0;l<3;l++) R1[j][k] += A0[j][l]*R0[l][k];
 899            } // end for j,k
 900        }else{
 901            for(j=0;j<3;j++) for(k=0;k<3;k++){
 902                R1[j][k] = 0.;
 903                for(l=0;l<3;l++) R1[j][k] += A1[j][l]*R0[l][k];
 904            } // end for j,k
 905        } // end if
 906        fGm[i]->SetMatrix(R1);
 907    } // end for i
 908    this->fTrans = (this->fTrans && 0xfffe) + 1;  // set bit 0 true.
 909    return;
 910 }
 911 //______________________________________________________________________
 912 Int_t AliITSgeom::GetNearest(const Double_t g[3],const Int_t lay=0){
 913 ////////////////////////////////////////////////////////////////////////
 914 //      Finds the Detector (Module) that is nearest the point g [cm] in
 915 // ALICE Global coordinates. If layer !=0 then the search is restricted
 916 // to Detectors (Modules) in that particular layer.
 917 ////////////////////////////////////////////////////////////////////////
 918      Int_t    i,l,a,e,in=0;
 919      Double_t d,dn=1.0e10;
 920      Bool_t   t=lay!=0; // skip if lay = 0 default value check all layers.
 921
 922      for(i=0;i<fNmodules;i++){
 923           if(t){GetModuleId(i,l,a,e);if(l!=lay) continue;}
 924           if((d=fGm[i]->Distance2(g))<dn){
 925                dn = d;
 926                in = i;
 927           } // end if
 928      } // end for i
 929      return in;
 930 }
 931 //______________________________________________________________________
 932 void AliITSgeom::GetNearest27(const Double_t g[3],Int_t n[27],const Int_t lay=0){
 933 ////////////////////////////////////////////////////////////////////////
 934 //      Finds 27 Detectors (Modules) that are nearest the point g [cm] in
 935 // ALICE Global coordinates. If layer !=0 then the search is restricted
 936 // to Detectors (Modules) in that particular layer. The number 27 comes
 937 // from including the nearest detector and all those around it (up, down,
 938 // left, right, forwards, backwards, and the corners).
 939 ////////////////////////////////////////////////////////////////////////
 940      Int_t    i,l,a,e,in[27]={0,0,0,0,0,0,0,0,0,
 941                               0,0,0,0,0,0,0,0,0,
 942                               0,0,0,0,0,0,0,0,0,};
 943      Double_t d,dn[27]={1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,
 944                         1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,
 945                         1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,
 946                         1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,1.0e10,
 947                         1.0e10,1.0e10,1.0e10};
 948      Bool_t   t=(lay!=0); // skip if lay = 0 default value check all layers.
 949
 950      for(i=0;i<fNmodules;i++){
 951           if(t){GetModuleId(i,l,a,e);if(l!=lay) continue;}
 952           for(a=0;a<27;a++){
 953                d = fGm[i]->Distance2(g);
 954                if(d<dn[a]){
 955                    for(e=26;e>a;e--){dn[e] = dn[e-1];in[e] = in[e-1];}
 956                    dn[a] = d; in[a] = i;
 957                } // end if d<dn[i]
 958           } // end for a
 959      } // end for i
 960      for(i=0;i<27;i++) n[i] = in[i];
 961 }
 962 //----------------------------------------------------------------------