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.3  1999/10/04 15:20:12  fca
  19 Correct syntax accepted by g++ but not standard for static members, remove minor warnings
  20
  21 Revision 1.2  1999/09/29 09:24:20  fca
  22 Introduction of the Copyright and cvs Log
  23
  24 */
  25
  26 ///////////////////////////////////////////////////////////////////////
  27 // ITS geometry manimulaiton routines.                               //
  28 // Created April 15 1999.                                            //
  29 // version: 0.0.0                                                    //
  30 // By: Bjorn S. Nilsen                                               //
  31 // version: 0.0.1                                                    //
  32 // Updated May 27 1999.                                              //
  33 // Added Cylinderical random and global based changes.               //
  34 // Added  function PrintComparison.                                  //
  35 ///////////////////////////////////////////////////////////////////////
  36 #include <iostream.h>
  37 #include <fstream.h>
  38 #include <iomanip.h>
  39 #include <stdio.h>
  40 #include "AliITSgeom.h"
  41 #include "TRandom.h"
  42
  43 ClassImp(AliITSgeom)
  44
  45 //_____________________________________________________________________
  46 AliITSgeom::AliITSgeom(){
  47 ////////////////////////////////////////////////////////////////////////
  48 //     The default constructor for the AliITSgeom class. It, by default,
  49 // sets fNlayers to zero and zeros all pointers.
  50 ////////////////////////////////////////////////////////////////////////
  51   // Default constructor.
  52   // Do not allocate anything zero everything
  53    fNlayers = 0;
  54    fNlad    = 0;
  55    fNdet    = 0;
  56    fg       = 0;
  57    fShape   = 0;
  58    return;
  59 }
  60
  61 //_____________________________________________________________________
  62 AliITSgeom::~AliITSgeom(){
  63 ////////////////////////////////////////////////////////////////////////
  64 //     The destructor for the AliITSgeom class. If the arrays fNlad,
  65 // fNdet, or fg have had memory allocated to them, there pointer values
  66 // are non zero, then this memory space is freed and they are set
  67 // to zero. In addition, fNlayers is set to zero. The destruction of
  68 // TObjArray fShape is, by default, handled by the TObjArray destructor.
  69 ////////////////////////////////////////////////////////////////////////
  70   // Default destructor.
  71   // if arrays exist delet them. Then set everything to zero.
  72    if(fg!=0){
  73       for(Int_t i=0;i<fNlayers;i++) delete[] fg[i];
  74       delete[] fg;
  75    } // end if fg!=0
  76    if(fNlad!=0) delete[] fNlad;
  77    if(fNdet!=0) delete[] fNdet;
  78    fNlayers = 0;
  79    fNlad    = 0;
  80    fNdet    = 0;
  81    fg       = 0;
  82    return;
  83 }
  84
  85 //_____________________________________________________________________
  86 AliITSgeom::AliITSgeom(const char *filename){
  87 ////////////////////////////////////////////////////////////////////////
  88 //     The constructor for the AliITSgeom class. All of the data to fill
  89 // this structure is read in from the file given my the input filename.
  90 ////////////////////////////////////////////////////////////////////////
  91    FILE     *pf;
  92    Int_t    i;
  93    ITS_geom *g;
  94    Int_t    l,a,d;
  95    Float_t  x,y,z,o,p,q,r,s,t;
  96    Double_t oor,pr,qr,rr,sr,tr; // Radians
  97    Double_t lr[9];
  98    Double_t si; // sin(angle)
  99    Double_t PI = TMath::Pi(), byPI = PI/180.;
 100
 101    pf = fopen(filename,"r");
 102
 103    fNlayers = 6; // set default number of ladders
 104    fNlad    = new Int_t[fNlayers];
 105    fNdet    = new Int_t[fNlayers];
 106    // find the number of laders and detectors in this geometry.
 107    for(i=0;i<fNlayers;i++){fNlad[i]=fNdet[i]=0;} // zero out arrays
 108    for(;;){ // for ever loop
 109       i = fscanf(pf,"%d %d %d %f %f %f %f %f %f %f %f %f",
 110                      &l,&a,&d,&x,&y,&z,&o,&p,&q,&r,&s,&t);
 111       if(i==EOF) break;
 112       if(l<1 || l>fNlayers) {
 113          printf("error in file %s layer=%d min is 1 max is %d/n",
 114                  filename,l,fNlayers);
 115          continue;
 116       }// end if l
 117       if(fNlad[l-1]<a) fNlad[l-1] = a;
 118       if(fNdet[l-1]<d) fNdet[l-1] = d;
 119    } // end for ever loop
 120    // counted the number of laders and detectors now allocate space.
 121    fg = new ITS_geom* [fNlayers];
 122    for(i=0;i<fNlayers;i++){
 123       fg[i] = 0;
 124       l = fNlad[i]*fNdet[i];
 125       fg[i] = new ITS_geom[l]; // allocate space for transforms
 126    } // end for i
 127
 128    // Set up Shapes for a default configuration of 6 layers.
 129    fShape = new TObjArray;
 130    AddShape((TObject *) new AliITSgeomSPD());  // shape 0
 131    AddShape((TObject *) new AliITSgeomSDD());  // shape 1
 132    AddShape((TObject *) new AliITSgeomSPD());  // shape 2
 133
 134    // prepair to read in transforms
 135    rewind(pf); // start over reading file
 136    for(;;){ // for ever loop
 137       i = fscanf(pf,"%d %d %d %f %f %f %f %f %f %f %f %f",
 138                      &l,&a,&d,&x,&y,&z,&o,&p,&q,&r,&s,&t);
 139       if(i==EOF) break;
 140       if(l<1 || l>fNlayers) {
 141          printf("error in file %s layer=%d min is 1 max is %d/n",
 142                  filename,l,fNlayers);
 143          continue;
 144       }// end if l
 145       l--; a--; d--; // shift layer, lader, and detector counters to zero base
 146       i = d + a*fNdet[l]; // position of this detector
 147       g = &(fg[l][i]);
 148
 149       oor = byPI*o;
 150       pr = byPI*p;
 151       qr = byPI*q;
 152       rr = byPI*r;
 153       sr = byPI*s;
 154       tr = byPI*t;
 155
 156       g->fx0   = x;
 157       g->fy0   = y;
 158       g->fz0   = z;
 159 //
 160       si    = sin(oor);if(o== 90.0) si = +1.0;
 161                       if(o==270.0) si = -1.0;
 162                       if(o==  0.0||o==180.) si = 0.0;
 163       lr[0] = si * cos(pr);
 164       lr[1] = si * sin(pr);
 165       lr[2] = cos(oor);if(o== 90.0||o==270.) lr[2] = 0.0;
 166                       if(o== 0.0)           lr[2] = +1.0;
 167                       if(o==180.0)          lr[2] = -1.0;
 168 //
 169       si    =  sin(qr);if(q== 90.0) si = +1.0;
 170                        if(q==270.0) si = -1.0;
 171                        if(q==  0.0||q==180.) si = 0.0;
 172       lr[3] = si * cos(rr);
 173       lr[4] = si * sin(rr);
 174       lr[5] = cos(qr);if(q== 90.0||q==270.) lr[5] = 0.0;
 175                       if(q==  0.0)          lr[5] = +1.0;
 176                       if(q==180.0)          lr[5] = -1.0;
 177 //
 178       si    = sin(sr);if(s== 90.0) si = +1.0;
 179                       if(s==270.0) si = -1.0;
 180                       if(s==  0.0||s==180.) si = 0.0;
 181       lr[6] = si * cos(tr);
 182       lr[7] = si * sin(tr);
 183       lr[8] = cos(sr);if(s== 90.0||s==270.0) lr[8] =  0.0;
 184                       if(s==  0.0)           lr[8] = +1.0;
 185                       if(s==180.0)           lr[8] = -1.0;
 186       // Normalize these elements
 187       for(a=0;a<3;a++){// reuse float si and integers a and d.
 188          si = 0.0;
 189          for(d=0;d<3;d++) si += lr[3*a+d]*lr[3*a+d];
 190          si = TMath::Sqrt(1./si);
 191          for(d=0;d<3;d++) g->fr[3*a+d] = lr[3*a+d] = si*lr[3*a+d];
 192       } // end for a
 193       // get angles from matrix up to a phase of 180 degrees.
 194       oor     = atan2(lr[7],lr[8]);if(oor<0.0) oor += 2.0*PI;
 195       pr     = asin(lr[2]);       if(pr<0.0) pr += 2.0*PI;
 196       qr     = atan2(lr[3],lr[0]);if(qr<0.0) qr += 2.0*PI;
 197       g->frx = oor;
 198       g->fry = pr;
 199       g->frz = qr;
 200       // l = layer-1 at this point.
 201            if(l==0||l==1) g->fShapeIndex = 0; // SPD's
 202       else if(l==2||l==3) g->fShapeIndex = 1; // SDD's
 203       else if(l==4||l==5) g->fShapeIndex = 2; // SSD's
 204    } // end for ever loop
 205    fclose(pf);
 206 }
 207
 208 //________________________________________________________________________
 209 AliITSgeom::AliITSgeom(AliITSgeom &source){
 210 ////////////////////////////////////////////////////////////////////////
 211 //     The copy constructor for the AliITSgeom class. It calls the
 212 // = operator function. See the = operator function for more details.
 213 ////////////////////////////////////////////////////////////////////////
 214    source = *this;  // Just use the = operator for now.
 215    return;
 216 }
 217
 218 //________________________________________________________________________
 219 void AliITSgeom::operator=(AliITSgeom &source){
 220 ////////////////////////////////////////////////////////////////////////
 221 //     The = operator function for the AliITSgeom class. It makes an
 222 // independent copy of the class in such a way that any changes made
 223 // to the copied class will not affect the source class in any way.
 224 // This is required for many ITS alignment studies where the copied
 225 // class is then modified by introducing some misalignment.
 226 ////////////////////////////////////////////////////////////////////////
 227    Int_t i,j,k;
 228
 229    if(this == &source) return; // don't assign to ones self.
 230
 231    // if there is an old structure allocated delete it first.
 232    if(fg != 0){
 233       for(i=0;i<fNlayers;i++) delete[] fg[i];
 234       delete[] fg;
 235    } // end if fg != 0
 236    if(fNlad != 0) delete[] fNlad;
 237    if(fNdet != 0) delete[] fNdet;
 238
 239    fNlayers = source.fNlayers;
 240    fNlad = new Int_t[fNlayers];
 241    for(i=0;i<fNlayers;i++) fNlad[i] = source.fNlad[i];
 242    fNdet = new Int_t[fNlayers];
 243    for(i=0;i<fNlayers;i++) fNdet[i] = source.fNdet[i];
 244    fShape = new TObjArray(*(source.fShape));//This does not make a proper copy.
 245    fg = new ITS_geom* [fNlayers];
 246    for(i=0;i<fNlayers;i++){
 247       fg[i] = new ITS_geom[fNlad[i]*fNdet[i]];
 248       for(j=0;j<(fNlad[i]*fNdet[i]);j++){
 249           fg[i][j].fShapeIndex = source.fg[i][j].fShapeIndex;
 250           fg[i][j].fx0 = source.fg[i][j].fx0;
 251           fg[i][j].fy0 = source.fg[i][j].fy0;
 252           fg[i][j].fz0 = source.fg[i][j].fz0;
 253           fg[i][j].frx = source.fg[i][j].frx;
 254           fg[i][j].fry = source.fg[i][j].fry;
 255           fg[i][j].frz = source.fg[i][j].frz;
 256           for(k=0;k<9;k++) fg[i][j].fr[k] = source.fg[i][j].fr[k];
 257       } // end for j
 258    } // end for i
 259    return;
 260 }
 261
 262
 263 //________________________________________________________________________
 264 void AliITSgeom::GtoL(Int_t lay,Int_t lad,Int_t det,
 265                        const Float_t *g,Float_t *l){
 266 ////////////////////////////////////////////////////////////////////////
 267 //     The function that does the global ALICE Cartesian coordinate
 268 // to local active volume detector Cartesian coordinate transformation.
 269 // The local detector coordinate system is determined by the layer,
 270 // ladder, and detector numbers. The global coordinates are entered by
 271 // the three element Float_t array g and the local coordinate values
 272 // are returned by the three element Float_t array l. The order of the
 273 // three elements are g[0]=x, g[1]=y, and g[2]=z, similarly for l.
 274 ////////////////////////////////////////////////////////////////////////
 275    Double_t x,y,z;
 276    ITS_geom *gl;
 277
 278    lay--; lad--; det--;
 279    gl = &(fg[lay][fNdet[lay]*lad+det]);
 280
 281    x    = g[0] - gl->fx0;
 282    y    = g[1] - gl->fy0;
 283    z    = g[2] - gl->fz0;
 284    l[0] = gl->fr[0]*x + gl->fr[1]*y + gl->fr[2]*z;
 285    l[1] = gl->fr[3]*x + gl->fr[4]*y + gl->fr[5]*z;
 286    l[2] = gl->fr[6]*x + gl->fr[7]*y + gl->fr[8]*z;
 287    return;
 288 }
 289
 290 //________________________________________________________________________
 291 void AliITSgeom::GtoL(const Int_t *id,const Float_t *g,Float_t *l){
 292 ////////////////////////////////////////////////////////////////////////
 293 //     The function that does the local active volume detector Cartesian
 294 // coordinate to global ALICE Cartesian coordinate transformation.
 295 // The local detector coordinate system is determined by the layer,
 296 // ladder, and detector numbers. The local coordinates are entered by
 297 // the three element Float_t array l and the global coordinate values
 298 // are returned by the three element Float_t array g. The order of the
 299 // three elements are l[0]=x, l[1]=y, and l[2]=z, similarly for g.
 300 ////////////////////////////////////////////////////////////////////////
 301    Int_t    lay,lad,det;
 302    Double_t x,y,z;
 303    ITS_geom *gl;
 304
 305    lay = id[0]; lad = id[1]; det = id[2];
 306    lay--; lad--; det--;
 307    gl = &(fg[lay][fNdet[lay]*lad+det]);
 308
 309    x    = g[0] - gl->fx0;
 310    y    = g[1] - gl->fy0;
 311    z    = g[2] - gl->fz0;
 312    l[0] = gl->fr[0]*x + gl->fr[1]*y + gl->fr[2]*z;
 313    l[1] = gl->fr[3]*x + gl->fr[4]*y + gl->fr[5]*z;
 314    l[2] = gl->fr[6]*x + gl->fr[7]*y + gl->fr[8]*z;
 315    return;
 316 }
 317 //________________________________________________________________________
 318 void AliITSgeom::GtoL(const Int_t index,const Float_t *g,Float_t *l){
 319 ////////////////////////////////////////////////////////////////////////
 320 //     The function that does the local active volume detector Cartesian
 321 // coordinate to global ALICE Cartesian coordinate transformation.
 322 // The local detector coordinate system is determined by the detector
 323 // index numbers (see GetModuleIndex and GetModuleID). The local
 324 // coordinates are entered by the three element Float_t array l and the
 325 // global coordinate values are returned by the three element Float_t array g.
 326 // The order of the three elements are l[0]=x, l[1]=y, and l[2]=z, similarly
 327 // for g.
 328 ////////////////////////////////////////////////////////////////////////
 329    Int_t    lay,lad,det;
 330    Double_t x,y,z;
 331    ITS_geom *gl;
 332
 333    this->GetModuleId(index,lay,lad,det);
 334    lay--; lad--; det--;
 335    gl = &(fg[lay][fNdet[lay]*lad+det]);
 336
 337    x    = g[0] - gl->fx0;
 338    y    = g[1] - gl->fy0;
 339    z    = g[2] - gl->fz0;
 340    l[0] = gl->fr[0]*x + gl->fr[1]*y + gl->fr[2]*z;
 341    l[1] = gl->fr[3]*x + gl->fr[4]*y + gl->fr[5]*z;
 342    l[2] = gl->fr[6]*x + gl->fr[7]*y + gl->fr[8]*z;
 343    return;
 344 }
 345
 346 //________________________________________________________________________
 347 void AliITSgeom::LtoG(Int_t lay,Int_t lad,Int_t det,
 348                       const Float_t *l,Float_t *g){
 349 ////////////////////////////////////////////////////////////////////////
 350 //     The function that does the local active volume detector Cartesian
 351 // coordinate to global ALICE Cartesian coordinate transformation.
 352 // The local detector coordinate system is determined by the layer,
 353 // ladder, and detector numbers. The local coordinates are entered by
 354 // the three element Float_t array l and the global coordinate values
 355 // are returned by the three element Float_t array g. The order of the
 356 // three elements are l[0]=x, l[1]=y, and l[2]=z, similarly for g.
 357 ////////////////////////////////////////////////////////////////////////
 358    Double_t x,y,z;
 359    ITS_geom *gl;
 360
 361    lay--; lad--; det--;
 362    gl   = &(fg[lay][fNdet[lay]*lad+det]);
 363
 364    x    = gl->fr[0]*l[0] + gl->fr[3]*l[1] + gl->fr[6]*l[2];
 365    y    = gl->fr[1]*l[0] + gl->fr[4]*l[1] + gl->fr[7]*l[2];
 366    z    = gl->fr[2]*l[0] + gl->fr[5]*l[1] + gl->fr[8]*l[2];
 367    g[0] = x + gl->fx0;
 368    g[1] = y + gl->fy0;
 369    g[2] = z + gl->fz0;
 370    return;
 371 }
 372
 373 //________________________________________________________________________
 374 void AliITSgeom::LtoG(const Int_t *id,const Float_t *l,Float_t *g){
 375 ////////////////////////////////////////////////////////////////////////
 376 //     The function that does the local active volume detector Cartesian
 377 // coordinate to global ALICE Cartesian coordinate transformation.
 378 // The local detector coordinate system is determined by the three
 379 // element array Id containing as it's three elements Id[0]=layer,
 380 // Id[1]=ladder, and Id[2]=detector numbers. The local coordinates
 381 // are entered by the three element Float_t array l and the global
 382 // coordinate values are returned by the three element Float_t array g.
 383 // The order of the three elements are l[0]=x, l[1]=y, and l[2]=z,
 384 // similarly for g.
 385 ////////////////////////////////////////////////////////////////////////
 386    Int_t    lay,lad,det;
 387    Double_t x,y,z;
 388    ITS_geom *gl;
 389
 390    lay = id[0]; lad = id[1]; det = id[2];
 391    lay--; lad--; det--;
 392    gl   = &(fg[lay][fNdet[lay]*lad+det]);
 393
 394    x    = gl->fr[0]*l[0] + gl->fr[3]*l[1] + gl->fr[6]*l[2];
 395    y    = gl->fr[1]*l[0] + gl->fr[4]*l[1] + gl->fr[7]*l[2];
 396    z    = gl->fr[2]*l[0] + gl->fr[5]*l[1] + gl->fr[8]*l[2];
 397    g[0] = x + gl->fx0;
 398    g[1] = y + gl->fy0;
 399    g[2] = z + gl->fz0;
 400    return;
 401 }
 402 //________________________________________________________________________
 403 void AliITSgeom::LtoG(const Int_t index,const Float_t *l,Float_t *g){
 404 ////////////////////////////////////////////////////////////////////////
 405 //     The function that does the local active volume detector Cartesian
 406 // coordinate to global ALICE Cartesian coordinate transformation.
 407 // The local detector coordinate system is determined by the detector
 408 // index number (see GetModuleIndex and GetModuleId). The local coordinates
 409 // are entered by the three element Float_t array l and the global
 410 // coordinate values are returned by the three element Float_t array g.
 411 // The order of the three elements are l[0]=x, l[1]=y, and l[2]=z,
 412 // similarly for g.
 413 ////////////////////////////////////////////////////////////////////////
 414    Int_t    lay,lad,det;
 415    Double_t x,y,z;
 416    ITS_geom *gl;
 417
 418    this->GetModuleId(index,lay,lad,det);
 419    lay--; lad--; det--;
 420    gl   = &(fg[lay][fNdet[lay]*lad+det]);
 421
 422    x    = gl->fr[0]*l[0] + gl->fr[3]*l[1] + gl->fr[6]*l[2];
 423    y    = gl->fr[1]*l[0] + gl->fr[4]*l[1] + gl->fr[7]*l[2];
 424    z    = gl->fr[2]*l[0] + gl->fr[5]*l[1] + gl->fr[8]*l[2];
 425    g[0] = x + gl->fx0;
 426    g[1] = y + gl->fy0;
 427    g[2] = z + gl->fz0;
 428    return;
 429 }
 430 //________________________________________________________________________
 431 void AliITSgeom::GtoLMomentum(Int_t lay,Int_t lad,Int_t det,
 432                               const Float_t *g,Float_t *l){
 433 ////////////////////////////////////////////////////////////////////////
 434 //     The function that does the global ALICE Cartesian momentum
 435 // to local active volume detector Cartesian momentum transformation.
 436 // The local detector coordinate system is determined by the layer,
 437 // ladder, and detector numbers. The global momentums are entered by
 438 // the three element Float_t array g and the local momentums values
 439 // are returned by the three element Float_t array l. The order of the
 440 // three elements are g[0]=x, g[1]=y, and g[2]=z, similarly for l.
 441 ////////////////////////////////////////////////////////////////////////
 442    Double_t px,py,pz;
 443    ITS_geom *gl;
 444
 445    lay--; lad--; det--;
 446    gl = &(fg[lay][fNdet[lay]*lad+det]);
 447
 448    px   = g[0];
 449    py   = g[1];
 450    pz   = g[2];
 451    l[0] = gl->fr[0]*px + gl->fr[1]*py + gl->fr[2]*pz;
 452    l[1] = gl->fr[3]*px + gl->fr[4]*py + gl->fr[5]*pz;
 453    l[2] = gl->fr[6]*px + gl->fr[7]*py + gl->fr[8]*pz;
 454    return;
 455 }
 456 //________________________________________________________________________
 457 void AliITSgeom::LtoGMomentum(Int_t lay,Int_t lad,Int_t det,
 458                               const Float_t *l,Float_t *g){
 459 ////////////////////////////////////////////////////////////////////////
 460 //     The function that does the local active volume detector Cartesian
 461 // momentum to global ALICE Cartesian momentum transformation.
 462 // The local detector momentum system is determined by the layer,
 463 // ladder, and detector numbers. The locall momentums are entered by
 464 // the three element Float_t array l and the global momentum values
 465 // are returned by the three element Float_t array g. The order of the
 466 // three elements are l[0]=x, l[1]=y, and l[2]=z, similarly for g.
 467 ////////////////////////////////////////////////////////////////////////
 468    Double_t px,py,pz;
 469    ITS_geom *gl;
 470
 471    lay--; lad--; det--;
 472    gl   = &(fg[lay][fNdet[lay]*lad+det]);
 473
 474    px   = gl->fr[0]*l[0] + gl->fr[3]*l[1] + gl->fr[6]*l[2];
 475    py   = gl->fr[1]*l[0] + gl->fr[4]*l[1] + gl->fr[7]*l[2];
 476    pz   = gl->fr[2]*l[0] + gl->fr[5]*l[1] + gl->fr[8]*l[2];
 477    g[0] = px;
 478    g[1] = py;
 479    g[2] = pz;
 480    return;
 481 }
 482 //___________________________________________________________________________
 483 Int_t AliITSgeom::GetModuleIndex(Int_t lay,Int_t lad,Int_t det){
 484     Int_t i,j,k;
 485
 486     i = fNdet[lay-1] * (lad-1) + det - 1;
 487     j = 0;
 488     for(k=0;k<lay-1;k++) j += fNdet[k]*fNlad[k];
 489     return (i+j);
 490 }
 491 //___________________________________________________________________________
 492 void AliITSgeom::GetModuleId(Int_t index,Int_t &lay,Int_t &lad,Int_t &det){
 493     Int_t i,j,k;
 494
 495     j = 0;
 496     for(k=0;k<fNlayers;k++){
 497         j += fNdet[k]*fNlad[k];
 498         if(j>index)break;
 499     } // end for k
 500     lay = k+1;
 501     i = index -j + fNdet[k]*fNlad[k];
 502     j = 0;
 503     for(k=0;k<fNlad[lay-1];k++){
 504         j += fNdet[lay-1];
 505         if(j>i)break;
 506     } // end for k
 507     lad = k+1;
 508     det = 1+i-fNdet[lay-1]*k;
 509     return;
 510 }
 511 //___________________________________________________________________________
 512 void AliITSgeom::GlobalChange(Float_t *tran,Float_t *rot){
 513 ////////////////////////////////////////////////////////////////////////
 514 //     This function performs a Cartesian translation and rotation of
 515 // the full ITS from its default position by an amount determined by
 516 // the three element arrays dtranslation and drotation. If every element
 517 // of dtranslation and drotation are zero then there is no change made
 518 // the geometry. The change is global in that the exact same translation
 519 // and rotation is done to every detector element in the exact same way.
 520 // The units of the translation are those of the Monte Carlo, usually cm,
 521 // and those of the rotation are in radians. The elements of dtranslation
 522 // are dtranslation[0] = x, dtranslation[1] = y, and dtranslation[2] = z.
 523 // The elements of drotation are drotation[0] = rx, drotation[1] = ry, and
 524 // drotation[2] = rz. A change in x will move the hole ITS in the ALICE
 525 // global x direction, the same for a change in y. A change in z will
 526 // result in a translation of the ITS as a hole up or down the beam line.
 527 // A change in the angles will result in the inclination of the ITS with
 528 // respect to the beam line, except for an effective rotation about the
 529 // beam axis which will just rotate the ITS as a hole about the beam axis.
 530 ////////////////////////////////////////////////////////////////////////
 531    Int_t    i,j,k,l;
 532    Double_t rx,ry,rz;
 533    Double_t sx,cx,sy,cy,sz,cz;
 534    ITS_geom *gl;
 535
 536    for(i=0;i<fNlayers;i++){
 537       for(j=0;j<fNlad[i];j++) for(k=0;k<fNdet[i];k++){
 538          l = fNdet[i]*j+k; // resolved index
 539          gl = &(fg[i][l]);
 540          gl->fx0 += tran[0];
 541          gl->fy0 += tran[1];
 542          gl->fz0 += tran[2];
 543          gl->frx +=  rot[0];
 544          gl->fry +=  rot[1];
 545          gl->frz +=  rot[2];
 546          rx = gl->frx; ry = gl->fry; rz = gl->frz;
 547          sx = sin(rx); cx = cos(rx);
 548          sy = sin(ry); cy = cos(ry);
 549          sz = sin(rz); cz = cos(rz);
 550          gl->fr[0] =  cz*cy;
 551          gl->fr[1] = -cz*sy*sx - sz*cx;
 552          gl->fr[2] = -cz*sy*cx + sz*sx;
 553          gl->fr[3] =  sz*cy;
 554          gl->fr[4] = -sz*sy*sx + cz*cx;
 555          gl->fr[5] = -sz*sy*cx - cz*sx;
 556          gl->fr[6] =  sy;
 557          gl->fr[7] =  cy*sx;
 558          gl->fr[8] =  cy*cx;
 559       } // end for j,k
 560    } // end for i
 561    return;
 562 }
 563
 564 //___________________________________________________________________________
 565 void AliITSgeom::GlobalCylindericalChange(Float_t *tran,Float_t *rot){
 566 ////////////////////////////////////////////////////////////////////////
 567 //     This function performs a cylindrical translation and rotation of
 568 // each ITS element by a fixed about in radius, rphi, and z from its
 569 // default position by an amount determined by the three element arrays
 570 // dtranslation and drotation. If every element of dtranslation and
 571 // drotation are zero then there is no change made the geometry. The
 572 // change is global in that the exact same distance change in translation
 573 // and rotation is done to every detector element in the exact same way.
 574 // The units of the translation are those of the Monte Carlo, usually cm,
 575 // and those of the rotation are in radians. The elements of dtranslation
 576 // are dtranslation[0] = r, dtranslation[1] = rphi, and dtranslation[2] = z.
 577 // The elements of drotation are drotation[0] = rx, drotation[1] = ry, and
 578 // drotation[2] = rz. A change in r will results in the increase of the
 579 // radius of each layer by the same about. A change in rphi will results in
 580 // the rotation of each layer by a different angle but by the same
 581 // circumferential distance. A change in z will result in a translation
 582 // of the ITS as a hole up or down the beam line. A change in the angles
 583 // will result in the inclination of the ITS with respect to the beam
 584 // line, except for an effective rotation about the beam axis which will
 585 // just rotate the ITS as a hole about the beam axis.
 586 ////////////////////////////////////////////////////////////////////////
 587    Int_t    i,j,k,l;
 588    Double_t rx,ry,rz,r,phi,rphi; // phi in radians
 589    Double_t sx,cx,sy,cy,sz,cz,r0;
 590    ITS_geom *gl;
 591
 592 //   printf("trans=%f %f %f rot=%f %f %f\n",tran[0],tran[1],tran[2],
 593 //        rot[0],rot[1],rot[2]);
 594    for(i=0;i<fNlayers;i++){
 595       for(j=0;j<fNlad[i];j++) for(k=0;k<fNdet[i];k++){
 596          l     = fNdet[i]*j+k; // resolved index
 597          gl    = &(fg[i][l]);
 598          r = r0= TMath::Hypot(gl->fy0,gl->fx0);
 599          phi   = atan2(gl->fy0,gl->fx0);
 600          rphi  = r0*phi;
 601          r    += tran[0];
 602          rphi += tran[1];
 603          phi   = rphi/r0;
 604          gl->fx0  = r*TMath::Cos(phi);
 605          gl->fy0  = r*TMath::Sin(phi);
 606          gl->fz0 += tran[2];
 607          gl->frx +=  rot[0];
 608          gl->fry +=  rot[1];
 609          gl->frz +=  rot[2];
 610          rx = gl->frx; ry = gl->fry; rz = gl->frz;
 611          sx = sin(rx); cx = cos(rx);
 612          sy = sin(ry); cy = cos(ry);
 613          sz = sin(rz); cz = cos(rz);
 614          gl->fr[0] =  cz*cy;
 615          gl->fr[1] = -cz*sy*sx - sz*cx;
 616          gl->fr[2] = -cz*sy*cx + sz*sx;
 617          gl->fr[3] =  sz*cy;
 618          gl->fr[4] = -sz*sy*sx + cz*cx;
 619          gl->fr[5] = -sz*sy*cx - cz*sx;
 620          gl->fr[6] =  sy;
 621          gl->fr[7] =  cy*sx;
 622          gl->fr[8] =  cy*cx;
 623       } // end for j,k
 624    } // end for i
 625    return;
 626 }
 627
 628 //___________________________________________________________________________
 629 void AliITSgeom::RandomChange(Float_t *stran,Float_t *srot){
 630 ////////////////////////////////////////////////////////////////////////
 631 //     This function performs a Gaussian random displacement and/or
 632 // rotation about the present global position of each active
 633 // volume/detector of the ITS. The sigma of the random displacement
 634 // is determined by the three element array stranslation, for the
 635 // x y and z translations, and the three element array srotation,
 636 // for the three rotation about the axis x y and z.
 637 ////////////////////////////////////////////////////////////////////////
 638    Int_t    i,j,k,l;
 639    Double_t rx,ry,rz;
 640    Double_t sx,cx,sy,cy,sz,cz;
 641    TRandom  ran;
 642    ITS_geom *gl;
 643
 644    for(i=0;i<fNlayers;i++){
 645       for(j=0;j<fNlad[i];j++) for(k=0;k<fNdet[i];k++){
 646          l = fNdet[i]*j+k; // resolved index
 647          gl = &(fg[i][l]);
 648          gl->fx0 += ran.Gaus(0.0,stran[0]);
 649          gl->fy0 += ran.Gaus(0.0,stran[1]);
 650          gl->fz0 += ran.Gaus(0.0,stran[2]);
 651          gl->frx += ran.Gaus(0.0, srot[0]);
 652          gl->fry += ran.Gaus(0.0, srot[1]);
 653          gl->frz += ran.Gaus(0.0, srot[2]);
 654          rx = gl->frx; ry = gl->fry; rz = gl->frz;
 655          sx = sin(rx); cx = cos(rx);
 656          sy = sin(ry); cy = cos(ry);
 657          sz = sin(rz); cz = cos(rz);
 658          gl->fr[0] =  cz*cy;
 659          gl->fr[1] = -cz*sy*sx - sz*cx;
 660          gl->fr[2] = -cz*sy*cx + sz*sx;
 661          gl->fr[3] =  sz*cy;
 662          gl->fr[4] = -sz*sy*sx + cz*cx;
 663          gl->fr[5] = -sz*sy*cx - cz*sx;
 664          gl->fr[6] =  sy;
 665          gl->fr[7] =  cy*sx;
 666          gl->fr[8] =  cy*cx;
 667       } // end for j,k
 668    } // end for i
 669    return;
 670 }
 671
 672 //___________________________________________________________________________
 673 void AliITSgeom::RandomCylindericalChange(Float_t *stran,Float_t *srot){
 674 ////////////////////////////////////////////////////////////////////////
 675 //     This function performs a Gaussian random displacement and/or
 676 // rotation about the present global position of each active
 677 // volume/detector of the ITS. The sigma of the random displacement
 678 // is determined by the three element array stranslation, for the
 679 // r rphi and z translations, and the three element array srotation,
 680 // for the three rotation about the axis x y and z. This random change
 681 // in detector position allow for the simulation of a random uncertainty
 682 // in the detector positions of the ITS.
 683 ////////////////////////////////////////////////////////////////////////
 684    Int_t     i,j,k,l;
 685    Double_t  rx,ry,rz,r,phi,x,y;  // phi in radians
 686    Double_t  sx,cx,sy,cy,sz,cz,r0;
 687    TRandom   ran;
 688    ITS_geom  *gl;
 689
 690 //   printf("trans=%f %f %f rot=%f %f %f\n",stran[0],stran[1],stran[2],
 691 //        srot[0],srot[1],srot[2]);
 692    for(i=0;i<fNlayers;i++){
 693       for(j=0;j<fNlad[i];j++) for(k=0;k<fNdet[i];k++){
 694          l     = fNdet[i]*j+k; // resolved index
 695          gl    = &(fg[i][l]);
 696          x     = gl->fx0;
 697          y     = gl->fy0;
 698          r = r0= TMath::Hypot(y,x);
 699          phi   = TMath::ATan2(y,x);
 700 //       if(phi<0.0) phi += 2.0*TMath::Pi();
 701          r    += ran.Gaus(0.0,stran[0]);
 702          phi  += ran.Gaus(0.0,stran[1])/r0;
 703 //       printf("fx0=%f fy0=%f rcos(phi)=%f rsin(phi)=%f\n",gl->fx0,gl->fy0,
 704 //              r*TMath::Cos(phi),r*TMath::Sin(phi));
 705          gl->fx0  = r*TMath::Cos(phi);
 706          gl->fy0  = r*TMath::Sin(phi);
 707 //       printf("r0=%f r=%f hypot=%f phi0=%f phi=%f ATan2=%f\n",
 708 //              r0,r,TMath::Hypot(gl->fy0,gl->fx0),
 709 //              phi0,phi,TMath::ATan2(gl->fy0,gl->fx0));
 710          gl->fz0 += ran.Gaus(0.0,stran[2]);
 711          gl->frx += ran.Gaus(0.0, srot[0]);
 712          gl->fry += ran.Gaus(0.0, srot[1]);
 713          gl->frz += ran.Gaus(0.0, srot[2]);
 714          rx = gl->frx; ry = gl->fry; rz = gl->frz;
 715          sx = sin(rx); cx = cos(rx);
 716          sy = sin(ry); cy = cos(ry);
 717          sz = sin(rz); cz = cos(rz);
 718          gl->fr[0] =  cz*cy;
 719          gl->fr[1] = -cz*sy*sx - sz*cx;
 720          gl->fr[2] = -cz*sy*cx + sz*sx;
 721          gl->fr[3] =  sz*cy;
 722          gl->fr[4] = -sz*sy*sx + cz*cx;
 723          gl->fr[5] = -sz*sy*cx - cz*sx;
 724          gl->fr[6] =  sy;
 725          gl->fr[7] =  cy*sx;
 726          gl->fr[8] =  cy*cx;
 727       } // end for j,k
 728    } // end for i
 729    return;
 730 }
 731
 732 //___________________________________________________________________________
 733 void AliITSgeom::SetByAngles(Int_t lay,Int_t lad,Int_t det,
 734                              Float_t rx,Float_t ry,Float_t rz){
 735 ////////////////////////////////////////////////////////////////////////
 736 //     This function computes a new rotation matrix based on the angles
 737 // rx, ry, and rz (in radians) for a give detector on the give ladder
 738 // in the give layer. A new
 739 // fg[layer-1][(fNlad[layer-1]*(ladder-1)+detector-1)].fr[] array is
 740 // computed.
 741 ////////////////////////////////////////////////////////////////////////
 742    ITS_geom *g;
 743    Double_t  sx,cx,sy,cy,sz,cz;
 744
 745    lay--; lad--; det--; // set to zero base now.
 746    g = &(fg[lay][fNdet[lay]*lad+det]);
 747
 748    sx = sin(rx); cx = cos(rx);
 749    sy = sin(ry); cy = cos(ry);
 750    sz = sin(rz); cz = cos(rz);
 751    g->frx   = rx;
 752    g->fry   = ry;
 753    g->frz   = rz;
 754    g->fr[0] =  cz*cy;
 755    g->fr[1] = -cz*sy*sx - sz*cx;
 756    g->fr[2] = -cz*sy*cx + sz*sx;
 757    g->fr[3] =  sz*cy;
 758    g->fr[4] = -sz*sy*sx + cz*cx;
 759    g->fr[5] = -sz*sy*cx - cz*sx;
 760    g->fr[6] =  sy;
 761    g->fr[7] =  cy*sx;
 762    g->fr[8] =  cy*cx;
 763    return;
 764 }
 765
 766 //___________________________________________________________________________
 767 void AliITSgeom::GetRotMatrix(Int_t lay,Int_t lad,Int_t det,Float_t *mat){
 768 ////////////////////////////////////////////////////////////////////////
 769 //     Returns, in the Float_t array pointed to by mat, the full rotation
 770 // matrix for the give detector defined by layer, ladder, and detector.
 771 // It returns all nine elements of fr in the ITS_geom structure. See the
 772 // description of the ITS_geom structure for further details of this
 773 // rotation matrix.
 774 ////////////////////////////////////////////////////////////////////////
 775    Int_t    i;
 776    ITS_geom *g;
 777
 778    lay--; lad--; det--; // shift to base 0
 779    g = &(fg[lay][fNdet[lay]*lad+det]);
 780    for(i=0;i<9;i++) mat[i] = g->fr[i];
 781    return;
 782 }
 783
 784 //___________________________________________________________________________
 785 void AliITSgeom::PrintComparison(FILE *fp,AliITSgeom *other){
 786 ////////////////////////////////////////////////////////////////////////
 787 //     This function was primarily created for diagnostic reasons. It
 788 // print to a file pointed to by the file pointer fp the difference
 789 // between two AliITSgeom classes. The format of the file is basicly,
 790 // define d? to be the difference between the same element of the two
 791 // classes. For example dfrx = this->fg[i][j].frx - other->fg[i][j].frx.
 792 // if(at least one of dfx0, dfy0, dfz0,dfrx,dfry,dfrz are non zero) then print
 793 // layer ladder detector dfx0 dfy0 dfz0 dfrx dfry dfrz
 794 // if(at least one of the 9 elements of dfr[] are non zero) then print
 795 // layer ladder detector dfr[0] dfr[1] dfr[2]
 796 //                       dfr[3] dfr[4] dfr[5]
 797 //                       dfr[6] dfr[7] dfr[8]
 798 // Only non zero values are printed to save space. The differences are
 799 // typical written to a file because there are usually a lot of numbers
 800 // printed out and it is usually easier to read them in some nice editor
 801 // rather than zooming quickly past you on a screen. fprintf is used to
 802 // do the printing. The fShapeIndex difference is not printed at this time.
 803 ////////////////////////////////////////////////////////////////////////
 804    Int_t    i,j,k,l;
 805    Double_t xt,yt,zt,xo,yo,zo;
 806    Double_t rxt,ryt,rzt,rxo,ryo,rzo;  // phi in radians
 807    ITS_geom *gt,*go;
 808    Bool_t   t;
 809
 810    for(i=0;i<this->fNlayers;i++){
 811       for(j=0;j<this->fNlad[i];j++) for(k=0;k<this->fNdet[i];k++){
 812          l   = this->fNdet[i]*j+k; // resolved index
 813          gt  = &(this->fg[i][l]);
 814          go  = &(other->fg[i][l]);
 815          xt  = gt->fx0; yt  = gt->fy0; zt  = gt->fz0;
 816          xo  = go->fx0; yo  = go->fy0; zo  = go->fz0;
 817          rxt = gt->frx; ryt = gt->fry; rzt = gt->frz;
 818          rxo = go->frx; ryo = go->fry; rzo = go->frz;
 819          if(!(xt==xo&&yt==yo&&zt==zo&&rxt==rxo&&ryt==ryo&&rzt==rzo))
 820          fprintf(fp,"%1.1d %2.2d %2.2d dTrans=%f %f %f drot=%f %f %f\n",
 821                  i+1,j+1,k+1,xt-xo,yt-yo,zt-zo,rxt-rxo,ryt-ryo,rzt-rzo);
 822          t = kFALSE;
 823          for(i=0;i<9;i++) t = gt->fr[i] != go->fr[i];
 824          if(t){
 825              fprintf(fp,"%1.1d %2.2d %2.2d dfr= %e %e %e\n",i+1,j+1,k+1,
 826                  gt->fr[0]-go->fr[0],gt->fr[1]-go->fr[1],gt->fr[2]-go->fr[2]);
 827              fprintf(fp,"        dfr= %e %e %e\n",
 828                  gt->fr[3]-go->fr[3],gt->fr[4]-go->fr[4],gt->fr[5]-go->fr[5]);
 829              fprintf(fp,"        dfr= %e %e %e\n",
 830                  gt->fr[6]-go->fr[6],gt->fr[7]-go->fr[7],gt->fr[8]-go->fr[8]);
 831          }
 832       } // end for j,k
 833    } // end for i
 834    return;
 835 }
 836
 837 //___________________________________________________________________________
 838 void AliITSgeom::PrintData(FILE *fp,Int_t lay,Int_t lad,Int_t det){
 839 ////////////////////////////////////////////////////////////////////////
 840 //     This function prints out the coordinate transformations for
 841 // the particular detector defined by layer, ladder, and detector
 842 // to the file pointed to by the File pointer fp. fprinf statements
 843 // are used to print out the numbers. The format is
 844 // layer ladder detector Trans= fx0 fy0 fz0 rot= frx fry frz Shape=fShapeIndex
 845 //                         dfr= fr[0] fr[1] fr[2]
 846 //                         dfr= fr[3] fr[4] fr[5]
 847 //                         dfr= fr[6] fr[7] fr[8]
 848 // By indicating which detector, some control over the information
 849 // is given to the user. The output it written to the file pointed
 850 // to by the file pointer fp. This can be set to stdout if you want.
 851 ////////////////////////////////////////////////////////////////////////
 852    Int_t    i,j,k,l;
 853    ITS_geom *gt;
 854
 855    i  = lay-1;
 856    j  = lad-1;
 857    k  = det-1;
 858    l  = this->fNdet[i]*j+k; // resolved index
 859    gt = &(this->fg[i][l]);
 860    fprintf(fp,"%1.1d %2.2d %2.2d Trans=%f %f %f rot=%f %f %f Shape=%d\n",
 861            i+1,j+1,k+1,gt->fx0,gt->fy0,gt->fz0,gt->frx,gt->fry,gt->frz,
 862            gt->fShapeIndex);
 863    fprintf(fp,"        dfr= %e %e %e\n",gt->fr[0],gt->fr[1],gt->fr[2]);
 864    fprintf(fp,"        dfr= %e %e %e\n",gt->fr[3],gt->fr[4],gt->fr[5]);
 865    fprintf(fp,"        dfr= %e %e %e\n",gt->fr[6],gt->fr[7],gt->fr[8]);
 866    return;
 867 }
 868 //___________________________________________________________________________
 869 void AliITSgeom::Streamer(TBuffer &R__b){
 870 ////////////////////////////////////////////////////////////////////////
 871 //     The default Streamer function "written by ROOT" doesn't write out
 872 // the arrays referenced by pointers. Therefore, a specific Streamer function
 873 // has to be written. This function should not be modified but instead added
 874 // on to so that older versions can still be read. The proper handling of
 875 // the version dependent streamer function hasn't been written do to the lack
 876 // of finding an example at the time of writting.
 877 ////////////////////////////////////////////////////////////////////////
 878    // Stream an object of class AliITSgeom.
 879     Int_t i,j,k;
 880
 881    if (R__b.IsReading()) {
 882       Version_t R__v = R__b.ReadVersion(); if (R__v) { }
 883       TObject::Streamer(R__b);
 884       R__b >> fNlayers;
 885       if(fNlad!=0) delete[] fNlad;
 886       if(fNdet!=0) delete[] fNdet;
 887       fNlad = new Int_t[fNlayers];
 888       fNdet = new Int_t[fNlayers];
 889       for(i=0;i<fNlayers;i++) R__b >> fNlad[i];
 890       for(i=0;i<fNlayers;i++) R__b >> fNdet[i];
 891       if(fg!=0){
 892           for(i=0;i<fNlayers;i++) delete[] fg[i];
 893           delete[] fg;
 894       } // end if fg!=0
 895       fg = new ITS_geom*[fNlayers];
 896       for(i=0;i<fNlayers;i++){
 897           fg[i] = new ITS_geom[fNlad[i]*fNdet[i]];
 898           for(j=0;j<fNlad[i]*fNdet[i];j++){
 899               R__b >> fg[i][j].fShapeIndex;
 900               R__b >> fg[i][j].fx0;
 901               R__b >> fg[i][j].fy0;
 902               R__b >> fg[i][j].fz0;
 903               R__b >> fg[i][j].frx;
 904               R__b >> fg[i][j].fry;
 905               R__b >> fg[i][j].frz;
 906               for(k=0;k<9;k++) R__b >> fg[i][j].fr[k];
 907           } // end for j
 908       } // end for i
 909       R__b >> fShape;
 910    } else {
 911       R__b.WriteVersion(AliITSgeom::IsA());
 912       TObject::Streamer(R__b);
 913       R__b << fNlayers;
 914       for(i=0;i<fNlayers;i++) R__b << fNlad[i];
 915       for(i=0;i<fNlayers;i++) R__b << fNdet[i];
 916       for(i=0;i<fNlayers;i++) for(j=0;j<fNlad[i]*fNdet[i];j++){
 917           R__b << fg[i][j].fShapeIndex;
 918           R__b << fg[i][j].fx0;
 919           R__b << fg[i][j].fy0;
 920           R__b << fg[i][j].fz0;
 921           R__b << fg[i][j].frx;
 922           R__b << fg[i][j].fry;
 923           R__b << fg[i][j].frz;
 924           for(k=0;k<9;k++) R__b << fg[i][j].fr[k];
 925       } // end for i,j
 926       R__b << fShape;
 927    }
 928 }
 929
 930 //___________________________________________________________________________
 931 ofstream & AliITSgeom::PrintGeom(ofstream &R__b){
 932 ////////////////////////////////////////////////////////////////////////
 933 //     The default Streamer function "written by ROOT" doesn't write out
 934 // the arrays referenced by pointers. Therefore, a specific Streamer function
 935 // has to be written. This function should not be modified but instead added
 936 // on to so that older versions can still be read. The proper handling of
 937 // the version dependent streamer function hasn't been written do to the lack
 938 // of finding an example at the time of writting.
 939 ////////////////////////////////////////////////////////////////////////
 940    // Stream an object of class AliITSgeom.
 941     Int_t i,j,k;
 942
 943     R__b.setf(ios::scientific);
 944     R__b << fNlayers << " ";
 945     for(i=0;i<fNlayers;i++) R__b << fNlad[i] << " ";
 946     for(i=0;i<fNlayers;i++) R__b << fNdet[i] << "\n";
 947     for(i=0;i<fNlayers;i++) for(j=0;j<fNlad[i]*fNdet[i];j++){
 948         R__b <<setprecision(16) << fg[i][j].fShapeIndex << " ";
 949         R__b <<setprecision(16) << fg[i][j].fx0 << " ";
 950         R__b <<setprecision(16) << fg[i][j].fy0 << " ";
 951         R__b <<setprecision(16) << fg[i][j].fz0 << " ";
 952         R__b <<setprecision(16) << fg[i][j].frx << " ";
 953         R__b <<setprecision(16) << fg[i][j].fry << " ";
 954         R__b <<setprecision(16) << fg[i][j].frz << "\n";
 955         for(k=0;k<9;k++) R__b <<setprecision(16) << fg[i][j].fr[k] << " ";
 956         R__b << "\n";
 957       } // end for i,j
 958 //      R__b << fShape;
 959       return R__b;
 960 }
 961
 962 //___________________________________________________________________________
 963 ifstream & AliITSgeom::ReadGeom(ifstream &R__b){
 964 ////////////////////////////////////////////////////////////////////////
 965 //     The default Streamer function "written by ROOT" doesn't write out
 966 // the arrays referenced by pointers. Therefore, a specific Streamer function
 967 // has to be written. This function should not be modified but instead added
 968 // on to so that older versions can still be read. The proper handling of
 969 // the version dependent streamer function hasn't been written do to the lack
 970 // of finding an example at the time of writting.
 971 ////////////////////////////////////////////////////////////////////////
 972    // Stream an object of class AliITSgeom.
 973     Int_t i,j,k;
 974
 975       R__b >> fNlayers;
 976       if(fNlad!=0) delete[] fNlad;
 977       if(fNdet!=0) delete[] fNdet;
 978       fNlad = new Int_t[fNlayers];
 979       fNdet = new Int_t[fNlayers];
 980       for(i=0;i<fNlayers;i++) R__b >> fNlad[i];
 981       for(i=0;i<fNlayers;i++) R__b >> fNdet[i];
 982       if(fg!=0){
 983           for(i=0;i<fNlayers;i++) delete[] fg[i];
 984           delete[] fg;
 985       } // end if fg!=0
 986       fg = new ITS_geom*[fNlayers];
 987       for(i=0;i<fNlayers;i++){
 988           fg[i] = new ITS_geom[fNlad[i]*fNdet[i]];
 989           for(j=0;j<fNlad[i]*fNdet[i];j++){
 990               R__b >> fg[i][j].fShapeIndex;
 991               R__b >> fg[i][j].fx0;
 992               R__b >> fg[i][j].fy0;
 993               R__b >> fg[i][j].fz0;
 994               R__b >> fg[i][j].frx;
 995               R__b >> fg[i][j].fry;
 996               R__b >> fg[i][j].frz;
 997               for(k=0;k<9;k++) R__b >> fg[i][j].fr[k];
 998           } // end for j
 999       } // end for i
1000 //      R__b >> fShape;
1001       return R__b;
1002 }