ITS/AliITSv11GeomCableRound.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
  18 // General Root includes
  19 //#include <Riostream.h>
  20 #include <TMath.h>
  21 #include <TVectorD.h>
  22
  23 // Root Geometry includes
  24 #include <TGeoManager.h>
  25 #include <TGeoVolume.h>
  26 #include <TGeoTube.h>
  27 #include <TGeoTorus.h>
  28 #include <TGeoMatrix.h>
  29
  30 #include "AliITSv11GeomCableRound.h"
  31
  32 //*************************************************************************
  33 //   Class for round cables
  34 //
  35 // Ludovic Gaudichet                                   gaudichet@to.infn.it
  36 //*************************************************************************
  37
  38 /*
  39 // ************************************************************************
  40 // Here is a example on how to use this class
  41 // ************************************************************************
  42
  43   // Getting some media
  44   TGeoMedium *air   = gGeoManager->GetMedium("ITS_AIR$");
  45   TGeoMedium *water = gGeoManager->GetMedium("ITS_WATER");
  46   TGeoMedium *alu   = gGeoManager->GetMedium("ITS_ITSal");
  47
  48   // Creating a small box inside a bigger one (containers)
  49   TGeoBBox *box1      = new TGeoBBox("box1", 6,10,10);
  50   TGeoBBox *bigBox    = new TGeoBBox("bigBox", 20,10,10);
  51   TGeoVolume *vbox1   = new TGeoVolume("vbox1", box1, air);
  52   TGeoVolume *vBigBox = new TGeoVolume("vBigBox", bigBox, air);
  53   vbox1->SetVisibility(kFALSE);
  54   vBigBox->SetVisibility(kFALSE);
  55
  56   TGeoTranslation *tr1 = new TGeoTranslation("negTr",-14,0,0);
  57   vBigBox->AddNode(vbox1, 1, tr1);
  58   moth->AddNode(vBigBox, 1, 0);
  59
  60   // **************************************************
  61   // Inserting a round cable (or here a water pipe...)
  62   // **************************************************
  63
  64   Int_t waterColor = 7;
  65   Int_t aluColor = 5;
  66   AliITSv11GeomCableRound roundCable("waterPipe", 0.9); //radius of 0.9cm
  67   roundCable.SetNLayers(2);
  68   roundCable.SetLayer(0, 0.7, water, waterColor); // radius of 0.7cm
  69   roundCable.SetLayer(1, 0.2, alu, aluColor);     // thickness of 0.2cm
  70
  71   // ****** Set check points and their containers ******
  72   // The 2 first points are in the small box (vbox1)
  73   // The second point is at the boundary
  74
  75   Double_t coord0[3] = {0,-2,-2};
  76   Double_t coord1[3] = {6,2,1};
  77   Double_t vect0[3]  = {1,1,0};
  78   Double_t vect1[3]  = {1,0,0};
  79   // coordinates have to be given in the specified container
  80   // reference system (here it's going to be vbox1).
  81   // vect1 and vect2 are vectors perpendicular to the segment ends
  82   // (These vectors don't need to be normalized)
  83   roundCable.AddCheckPoint( vbox1, 0, coord0, vect0);
  84   roundCable.AddCheckPoint( vbox1, 1, coord1, vect1);
  85
  86   // Then, let's cross the boundary ! You just need
  87   // to put the next point in the other volume, vBigBox.
  88   // At the moment of creating the second segment, it will
  89   // be inserted in this volume. That is why the point 1 had to
  90   // be at the boundary, because otherwise the second segment
  91   // between de points 1 and 2 would have been inserted in the
  92   // vBigBox but in the same time would have cross its
  93   // boundary ...
  94   Double_t coord2[3] = {-2,6,4}; // coord. syst. of vBigBox !
  95   Double_t vect2[3]= {1,1,0.5};
  96   roundCable.AddCheckPoint( vBigBox, 2, coord2, vect2);
  97
  98   Double_t coord3[3] = {4,6,4};
  99   Double_t vect3[3]= {-1,0,0};
 100   roundCable.AddCheckPoint( vBigBox, 3, coord3, vect3);
 101
 102   Double_t coord4[3] = {4,0,-4};
 103   Double_t vect4[3]= {1,0,0};
 104   roundCable.AddCheckPoint( vBigBox, 4, coord4, vect4);
 105
 106   Double_t coord5[3] = {4,-6,4};
 107   Double_t vect5[3]= {1,0,0};
 108   roundCable.AddCheckPoint( vBigBox, 5, coord5, vect5);
 109
 110   Double_t coord6[3] = {7,-6,4};
 111   Double_t vect6[3]= {1,0,0};
 112   roundCable.AddCheckPoint( vBigBox, 6, coord6, vect6);
 113
 114   Double_t r = 7;
 115   Double_t angle = 70*TMath::DegToRad();
 116   Double_t coord7[3] = {coord6[0] +r*sin(angle), coord6[1],
 117                         coord6[2] -r*(1-cos(angle)) };
 118   Double_t vect7[3]= {r*cos(angle),0,-r*sin(angle)};
 119   roundCable.AddCheckPoint( vBigBox, 7, coord7, vect7);
 120
 121   Double_t coord8[3] = { coord7[0]+vect7[0], coord7[1]+vect7[1],-10};
 122   Double_t vect8[3]= {0,0,1};
 123   roundCable.AddCheckPoint( vBigBox, 8, coord8, vect8);
 124
 125   // ****** Creating the corresponding volume ******
 126   // Since the container volumes of the check points have
 127   // been recorded, this can be done at any moments, providing
 128   // that the container volumes are found in the sub-nodes
 129   // of the initial node (the top volume of the TGeoManager or
 130   // the volume set in SetInitialNode(TGeoVolume*) function)
 131
 132   roundCable.SetInitialNode(vBigBox); //Set the root node
 133   roundCable.CreateAndInsertCableSegment( 1);
 134   // This command means : create the segment between point 0
 135   // and point 1. The segment is automatically inserted in the
 136   // container volume of point 1.
 137
 138   roundCable.CreateAndInsertCableSegment( 2);
 139   roundCable.CreateAndInsertCableSegment( 3);
 140
 141   // The following segment is going to be a torus segment.
 142   // The radius and position of the torus is defined by the
 143   // orthogonal vector of point 4 (the orientation of this vector
 144   // and the position of the 2 check points are enough to define
 145   // completely the torus)
 146   roundCable.CreateAndInsertTorusSegment( 4, 180);
 147   // The second argument is an additionnal rotation of the
 148   // segment around the axis defined by the 2 check points.
 149
 150   roundCable.CreateAndInsertTorusSegment( 5);
 151   roundCable.CreateAndInsertCableSegment( 6);
 152   roundCable.CreateAndInsertTorusSegment( 7,180);
 153   roundCable.CreateAndInsertCableSegment( 8);
 154
 155 */
 156
 157
 158
 159 ClassImp(AliITSv11GeomCableRound)
 160
 161 //________________________________________________________________________
 162 AliITSv11GeomCableRound::
 163 AliITSv11GeomCableRound(const char* name, Double_t radius) :
 164   AliITSv11GeomCable(name),
 165   fRadius(radius),
 166   fNlayer(0),
 167   fPhiMin(0),
 168   fPhiMax(360)
 169  {
 170   // Constructor
 171   for (Int_t i=0; i<fgkCableMaxLayer ; i++) {
 172     fLayThickness[i] = 0;
 173     fLayColor[i] = 0;
 174     fLayMedia[i] = 0;
 175   };
 176 }
 177
 178 //________________________________________________________________________
 179 AliITSv11GeomCableRound::AliITSv11GeomCableRound(const AliITSv11GeomCableRound &s) :
 180   AliITSv11GeomCable(s),fRadius(s.fRadius),fNlayer(s.fNlayer),fPhiMin(s.fPhiMin),
 181   fPhiMax(s.fPhiMax)
 182 {
 183   //     Copy Constructor
 184   for (Int_t i=0; i<s.fNlayer; i++) {
 185     fLayThickness[i] = s.fLayThickness[i];
 186     fLayMedia[i] = s.fLayMedia[i];
 187     fLayColor[i] = s.fLayColor[i];
 188   }
 189 }
 190
 191 //________________________________________________________________________
 192 AliITSv11GeomCableRound& AliITSv11GeomCableRound::
 193 operator=(const AliITSv11GeomCableRound &s) {
 194   //     Assignment operator
 195   // Not fully inplemented yet !!!
 196
 197   if(&s == this) return *this;
 198   *this = s;
 199   fRadius = s.fRadius;
 200   fPhiMin = s.fPhiMin;
 201   fPhiMax = s.fPhiMax;
 202   fNlayer = s.fNlayer;
 203   for (Int_t i=0; i<s.fNlayer; i++) {
 204     fLayThickness[i] = s.fLayThickness[i];
 205     fLayMedia[i] = s.fLayMedia[i];
 206     fLayColor[i] = s.fLayColor[i];
 207   };
 208   return *this;
 209 }
 210
 211 //________________________________________________________________________
 212 Int_t AliITSv11GeomCableRound::GetPoint( Int_t iCheckPt, Double_t *coord)
 213   const {
 214   // Get check point #iCheckPt
 215   TVectorD *coordVector =(TVectorD *)fPointArray.UncheckedAt(2*iCheckPt);
 216 #if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
 217   CopyFrom(coord, coordVector->GetElements());
 218 #else
 219   CopyFrom(coord, coordVector->GetMatrixArray());
 220 #endif
 221   return kTRUE;
 222 }
 223
 224 //________________________________________________________________________
 225 Int_t AliITSv11GeomCableRound::GetVect( Int_t iCheckPt, Double_t *coord)
 226   const {
 227   //
 228   // Get vector transverse to the section at point #iCheckPt
 229   //
 230
 231   TVectorD *coordVector =(TVectorD *)fPointArray.UncheckedAt(2*iCheckPt+1);
 232 #if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
 233   CopyFrom(coord, coordVector->GetElements());
 234 #else
 235   CopyFrom(coord, coordVector->GetMatrixArray());
 236 #endif
 237   return kTRUE;
 238 }
 239
 240 //________________________________________________________________________
 241 void AliITSv11GeomCableRound::AddCheckPoint( TGeoVolume *vol, Int_t iCheckPt,
 242                                                Double_t *coord, Double_t *orthVect)
 243 {
 244   //
 245   // Add point #iCheckPt and its transverse vector. Point is added at (i) in
 246   // fPointArray and the vector is added at (i+1)
 247   //
 248
 249
 250   if (iCheckPt>=fVolumeArray.GetEntriesFast()) {
 251     fVolumeArray.AddLast(vol);
 252     TVectorD *point = new TVectorD(3,coord);
 253     TVectorD *vect  = new TVectorD(3,orthVect);
 254     fPointArray.AddLast(point);
 255     fPointArray.AddLast(vect);
 256
 257   } else if ((iCheckPt >= 0)&&(iCheckPt < fVolumeArray.GetEntriesFast())) {
 258     fVolumeArray.AddAt(vol, iCheckPt);
 259     TVectorD *point = new TVectorD(3,coord);
 260     TVectorD *vect  = new TVectorD(3,orthVect);
 261     fPointArray.AddAt(point, iCheckPt*2  );
 262     fPointArray.AddAt(vect,  iCheckPt*2+1);
 263   };
 264 }
 265
 266 //________________________________________________________________________
 267 void AliITSv11GeomCableRound::PrintCheckPoints() const {
 268   // Print all check points
 269
 270   printf("  ---\n  Printing all check points of the round cable\n");
 271   for (Int_t i = 0; i<fVolumeArray.GetEntriesFast(); i++) {
 272     TVectorD *coordVector = (TVectorD *)fPointArray.UncheckedAt(i*2);
 273     //TVectorD *vectVector = (TVectorD *)fPointArray.UncheckedAt(i*2+1);
 274     Double_t coord[3];
 275 #if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
 276     CopyFrom(coord, coordVector->GetElements());
 277 #else
 278     CopyFrom(coord, coordVector->GetMatrixArray());
 279 #endif
 280     printf("   ( %.2f, %.2f, %.2f )\n", coord[0], coord[1], coord[2]);
 281   };
 282 }
 283
 284
 285 //________________________________________________________________________
 286 TGeoVolume* AliITSv11GeomCableRound::CreateAndInsertCableSegment(Int_t p2,
 287                                                                  TGeoCombiTrans** ct)
 288 {
 289 //    Creates a cable segment between points p1 and p2.
 290 //
 291 // The segment volume is created inside the volume containing point2
 292 // Therefore this segment should be defined in this volume only.
 293 // I mean here that, if the previous point is in another volume,
 294 // it should be just at the border between the 2 volumes. Also the
 295 // orientation vector of the previous point should be othogonal to
 296 // the surface between the 2 volumes.
 297
 298   TGeoNode *mainNode;
 299   if (fInitialNode==0) {
 300     TObjArray *nodes = gGeoManager->GetListOfNodes();
 301     if (nodes->GetEntriesFast()==0) return 0;
 302     mainNode = (TGeoNode *) nodes->UncheckedAt(0);
 303   } else {
 304     mainNode = fInitialNode;
 305   };
 306
 307   Int_t p1 = p2 - 1;
 308   TGeoVolume *p1Vol = GetVolume(p1);
 309   TGeoVolume *p2Vol = GetVolume(p2);
 310
 311   ResetCheckDaughter();
 312   fCurrentVol = p1Vol;
 313   if (! CheckDaughter(mainNode)) {
 314     printf("Error::volume containing point is not visible in node tree!\n");
 315     return 0;
 316   };
 317
 318   Double_t coord1[3], coord2[3], vect1[3], vect2[3];
 319   //=================================================
 320   // Get p1 position in the systeme of p2
 321   if (p1Vol!=p2Vol) {
 322
 323     Int_t p1nodeInd[fgkCableMaxNodeLevel];
 324     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
 325     Int_t p1volLevel = 0;
 326     while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
 327     p1volLevel--;
 328
 329     ResetCheckDaughter();
 330     fCurrentVol = p2Vol;
 331     if (! CheckDaughter(mainNode)) {
 332       printf("Error::volume containing point is not visible in node tree!\n");
 333       return 0;
 334     };
 335     Int_t p2nodeInd[fgkCableMaxNodeLevel];
 336     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
 337     Int_t commonMotherLevel = 0;
 338     while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
 339       commonMotherLevel++;
 340     commonMotherLevel--;
 341     Int_t p2volLevel = 0;
 342     while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
 343     p2volLevel--;
 344
 345     // Get coord and vect of p1 in the common mother reference system
 346     GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1);
 347     GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1);
 348     // Translate them in the reference system of the volume containing p2
 349     TGeoNode *pathNode[fgkCableMaxNodeLevel];
 350     pathNode[0] = mainNode;
 351     for (Int_t i=0; i<=p2volLevel; i++) {
 352       pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
 353     };
 354     Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
 355     Double_t globalVect1[3]  = {vect1[0], vect1[1], vect1[2]};
 356
 357     for (Int_t i = commonMotherLevel+1; i<=p2volLevel; i++) {
 358       pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
 359       pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
 360       CopyFrom(globalCoord1, coord1);
 361       CopyFrom(globalVect1, vect1);
 362     };
 363   } else {
 364     GetCheckPoint(p1, 0, 0, coord1);
 365     GetCheckVect(p1, 0, 0, vect1);
 366   };
 367
 368   //=================================================
 369   // Get p2 position in the systeme of p2
 370   GetCheckPoint(p2, 0, 0, coord2);
 371   GetCheckVect(p2, 0, 0, vect2);
 372
 373   Double_t cx = (coord1[0]+coord2[0])/2;
 374   Double_t cy = (coord1[1]+coord2[1])/2;
 375   Double_t cz = (coord1[2]+coord2[2])/2;
 376   Double_t dx = coord2[0]-coord1[0];
 377   Double_t dy = coord2[1]-coord1[1];
 378   Double_t dz = coord2[2]-coord1[2];
 379
 380   //=================================================
 381   // Positionning of the segment between the 2 points
 382   if ((dy<1e-31)&&(dy>0)) dy = 1e-31;
 383   if ((dz<1e-31)&&(dz>0)) dz = 1e-31;
 384   if ((dy>-1e-31)&&(dy<0)) dy = -1e-31;
 385   if ((dz>-1e-31)&&(dz<0)) dz = -1e-31;
 386
 387   Double_t angleRot1 = -TMath::ATan2(dx,dy);
 388   Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
 389   Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
 390   TGeoRotation *rot = new TGeoRotation("",angleRot1*TMath::RadToDeg(),
 391                                        angleRotDiag*TMath::RadToDeg(),
 392                                        0);
 393   Double_t localVect1[3], localVect2[3];
 394   rot->MasterToLocalVect(vect1, localVect1);
 395   rot->MasterToLocalVect(vect2, localVect2);
 396   TGeoTranslation *trans = new TGeoTranslation("",cx, cy, cz);
 397
 398   //=================================================
 399   // Create the segment and add it to the mother volume
 400   TGeoVolume *vCableSeg = CreateSegment(coord1, coord2,
 401                                         localVect1, localVect2, p2);
 402
 403   TGeoCombiTrans  *combi = new TGeoCombiTrans(*trans, *rot);
 404   p2Vol->AddNode(vCableSeg, p2, combi);
 405   //=================================================
 406   delete rot;
 407   delete trans;
 408
 409   if (fDebug) {
 410     printf("---\n  Cable segment points : ");
 411     printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
 412     printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
 413   };
 414 //   #include <TGeoSphere.h>
 415 //   TGeoMedium *airSDD = gGeoManager->GetMedium("ITS_AIR$");
 416 //   TGeoSphere *sphere = new TGeoSphere(0, 0.15);
 417 //   TGeoVolume *vSphere = new TGeoVolume("", sphere, airSDD);
 418 //   TGeoTranslation *trC = new TGeoTranslation("", cx, cy, cz);
 419 //   TGeoTranslation *tr1 = new TGeoTranslation("",coord1[0],
 420 //                                           coord1[1],coord1[2]);
 421 //   TGeoTranslation *tr2 = new TGeoTranslation("",coord2[0],
 422 //                                           coord2[1],coord2[2]);
 423 //   p2Vol->AddNode(vSphere, p2*3-2, trC);
 424 //   p2Vol->AddNode(vSphere, p2*3-1, tr1);
 425 //   p2Vol->AddNode(vSphere, p2*3  , tr2);
 426
 427   if (ct) *ct = combi;
 428   return vCableSeg;
 429 }
 430
 431 //________________________________________________________________________
 432 TGeoVolume* AliITSv11GeomCableRound::CreateAndInsertTubeSegment(Int_t p2,
 433                                                                  TGeoCombiTrans** ct)
 434 {
 435 //    Creates a cable segment between points p1 and p2.
 436 //
 437 //  This creates simple tube sections, i.e. the cable ends are
 438 // cutted perpendicularly to the tube axis. The method has to
 439 // be used only in this simple case, in ordder to save some memory
 440
 441   TGeoNode *mainNode;
 442   if (fInitialNode==0) {
 443     TObjArray *nodes = gGeoManager->GetListOfNodes();
 444     if (nodes->GetEntriesFast()==0) return 0;
 445     mainNode = (TGeoNode *) nodes->UncheckedAt(0);
 446   } else {
 447     mainNode = fInitialNode;
 448   };
 449
 450   Int_t p1 = p2 - 1;
 451   TGeoVolume *p1Vol = GetVolume(p1);
 452   TGeoVolume *p2Vol = GetVolume(p2);
 453
 454   ResetCheckDaughter();
 455   fCurrentVol = p1Vol;
 456   if (! CheckDaughter(mainNode)) {
 457     printf("Error::volume containing point is not visible in node tree!\n");
 458     return 0;
 459   };
 460
 461   Double_t coord1[3], coord2[3], vect1[3], vect2[3];
 462   //=================================================
 463   // Get p1 position in the systeme of p2
 464   if (p1Vol!=p2Vol) {
 465
 466     Int_t p1nodeInd[fgkCableMaxNodeLevel];
 467     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
 468     Int_t p1volLevel = 0;
 469     while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
 470     p1volLevel--;
 471
 472     ResetCheckDaughter();
 473     fCurrentVol = p2Vol;
 474     if (! CheckDaughter(mainNode)) {
 475       printf("Error::volume containing point is not visible in node tree!\n");
 476       return 0;
 477     };
 478     Int_t p2nodeInd[fgkCableMaxNodeLevel];
 479     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
 480     Int_t commonMotherLevel = 0;
 481     while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
 482       commonMotherLevel++;
 483     commonMotherLevel--;
 484     Int_t p2volLevel = 0;
 485     while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
 486     p2volLevel--;
 487
 488     // Get coord and vect of p1 in the common mother reference system
 489     GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1);
 490     GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1);
 491     // Translate them in the reference system of the volume containing p2
 492     TGeoNode *pathNode[fgkCableMaxNodeLevel];
 493     pathNode[0] = mainNode;
 494     for (Int_t i=0; i<=p2volLevel; i++) {
 495       pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
 496     };
 497     Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
 498     Double_t globalVect1[3]  = {vect1[0], vect1[1], vect1[2]};
 499
 500     for (Int_t i = commonMotherLevel+1; i<=p2volLevel; i++) {
 501       pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
 502       pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
 503       CopyFrom(globalCoord1, coord1);
 504       CopyFrom(globalVect1, vect1);
 505     };
 506   } else {
 507     GetCheckPoint(p1, 0, 0, coord1);
 508     GetCheckVect(p1, 0, 0, vect1);
 509   };
 510
 511   //=================================================
 512   // Get p2 position in the systeme of p2
 513   GetCheckPoint(p2, 0, 0, coord2);
 514   GetCheckVect(p2, 0, 0, vect2);
 515
 516   Double_t cx = (coord1[0]+coord2[0])/2;
 517   Double_t cy = (coord1[1]+coord2[1])/2;
 518   Double_t cz = (coord1[2]+coord2[2])/2;
 519   Double_t dx = coord2[0]-coord1[0];
 520   Double_t dy = coord2[1]-coord1[1];
 521   Double_t dz = coord2[2]-coord1[2];
 522
 523   //=================================================
 524   // Positionning of the segment between the 2 points
 525   if ((dy<1e-31)&&(dy>0)) dy = 1e-31;
 526   if ((dz<1e-31)&&(dz>0)) dz = 1e-31;
 527   if ((dy>-1e-31)&&(dy<0)) dy = -1e-31;
 528   if ((dz>-1e-31)&&(dz<0)) dz = -1e-31;
 529
 530   Double_t angleRot1 = -TMath::ATan2(dx,dy);
 531   Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
 532   Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
 533   TGeoRotation *rot = new TGeoRotation("",angleRot1*TMath::RadToDeg(),
 534                                        angleRotDiag*TMath::RadToDeg(),
 535                                        0);
 536   TGeoTranslation *trans = new TGeoTranslation("",cx, cy, cz);
 537
 538   //=================================================
 539   // Create the segment and add it to the mother volume
 540   TGeoVolume *vCableSeg = CreateTubeSegment( coord1,coord2, p2);
 541
 542   TGeoCombiTrans  *combi = new TGeoCombiTrans(*trans, *rot);
 543   p2Vol->AddNode(vCableSeg, p2, combi);
 544   //=================================================
 545   delete rot;
 546   delete trans;
 547
 548   if (fDebug) {
 549     printf("---\n  Cable segment points : ");
 550     printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
 551     printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
 552   };
 553
 554   if (ct) *ct = combi;
 555   return vCableSeg;
 556 }
 557
 558 //________________________________________________________________________
 559 TGeoVolume* AliITSv11GeomCableRound::CreateAndInsertTorusSegment(Int_t p2,
 560                                                                  Double_t rotation,
 561                                                                  TGeoCombiTrans** ct)
 562 {
 563   // Create a torus cable segment between points p1 and p2.
 564   // The radius and position of the torus is defined by the
 565   // perpendicular vector of point p2 (the orientation of this vector
 566   // and the position of the 2 check points are enough to completely
 567   // define the torus)
 568
 569   TGeoNode *mainNode;
 570   if (fInitialNode==0) {
 571     TObjArray *nodes = gGeoManager->GetListOfNodes();
 572     if (nodes->GetEntriesFast()==0) return 0;
 573     mainNode = (TGeoNode *) nodes->UncheckedAt(0);
 574   } else {
 575     mainNode = fInitialNode;
 576   };
 577
 578   Int_t p1 = p2 - 1;
 579   TGeoVolume *p1Vol = GetVolume(p1);
 580   TGeoVolume *p2Vol = GetVolume(p2);
 581
 582   ResetCheckDaughter();
 583   fCurrentVol = p1Vol;
 584   if (! CheckDaughter(mainNode)) {
 585     printf("Error::volume containing point is not visible in node tree!\n");
 586     return 0;
 587   };
 588
 589   Double_t coord1[3], coord2[3], vect1[3], vect2[3];
 590   //=================================================
 591   // Get p1 position in the systeme of p2
 592   if (p1Vol!=p2Vol) {
 593
 594     Int_t p1nodeInd[fgkCableMaxNodeLevel];
 595     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
 596     Int_t p1volLevel = 0;
 597     while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
 598     p1volLevel--;
 599
 600     ResetCheckDaughter();
 601     fCurrentVol = p2Vol;
 602     if (! CheckDaughter(mainNode)) {
 603       printf("Error::volume containing point is not visible in node tree!\n");
 604       return 0;
 605     };
 606     Int_t p2nodeInd[fgkCableMaxNodeLevel];
 607     for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
 608     Int_t commonMotherLevel = 0;
 609     while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
 610       commonMotherLevel++;
 611     commonMotherLevel--;
 612     Int_t p2volLevel = 0;
 613     while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
 614     p2volLevel--;
 615
 616     // Get coord and vect of p1 in the common mother reference system
 617     GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1);
 618     GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1);
 619     // Translate them in the reference system of the volume containing p2
 620     TGeoNode *pathNode[fgkCableMaxNodeLevel];
 621     pathNode[0] = mainNode;
 622     for (Int_t i=0; i<=p2volLevel; i++) {
 623       pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
 624     };
 625     Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
 626     Double_t globalVect1[3]  = {vect1[0], vect1[1], vect1[2]};
 627
 628     for (Int_t i = commonMotherLevel+1; i<=p2volLevel; i++) {
 629       pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
 630       pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
 631       CopyFrom(globalCoord1, coord1);
 632       CopyFrom(globalVect1, vect1);
 633     };
 634   } else {
 635     GetCheckPoint(p1, 0, 0, coord1);
 636     GetCheckVect(p1, 0, 0, vect1);
 637   };
 638
 639   //=================================================
 640   // Get p2 position in the systeme of p2
 641   GetCheckPoint(p2, 0, 0, coord2);
 642   GetCheckVect(p2, 0, 0, vect2);
 643
 644   Double_t cx = (coord1[0]+coord2[0])/2;
 645   Double_t cy = (coord1[1]+coord2[1])/2;
 646   Double_t cz = (coord1[2]+coord2[2])/2;
 647   Double_t dx = coord2[0]-coord1[0];
 648   Double_t dy = coord2[1]-coord1[1];
 649   Double_t dz = coord2[2]-coord1[2];
 650   Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
 651
 652   //=================================================
 653   // Positionning of the segment between the 2 points
 654   if ((dy<1e-31)&&(dy>0)) dy = 1e-31;
 655   if ((dz<1e-31)&&(dz>0)) dz = 1e-31;
 656   if ((dy>-1e-31)&&(dy<0)) dy = -1e-31;
 657   if ((dz>-1e-31)&&(dz<0)) dz = -1e-31;
 658
 659   Double_t angleRot1 = -TMath::ATan2(dx,dy);
 660   Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
 661   Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
 662
 663   TGeoRotation rotTorusTemp("",angleRot1*TMath::RadToDeg(),
 664                             angleRotDiag*TMath::RadToDeg(),0);
 665   TGeoRotation rotTorusToZ("",0,90,0);
 666   rotTorusTemp.MultiplyBy(&rotTorusToZ, kTRUE);
 667   Double_t localVect2[3];
 668   rotTorusTemp.MasterToLocalVect(vect2, localVect2);
 669   if (localVect2[1]<0) {
 670     localVect2[0] = -localVect2[0];
 671     localVect2[1] = -localVect2[1];
 672     localVect2[2] = -localVect2[2];
 673   };
 674   Double_t normVect2 = TMath::Sqrt(localVect2[0]*localVect2[0]+
 675                                    localVect2[1]*localVect2[1]+
 676                                    localVect2[2]*localVect2[2]);
 677   Double_t axisX[3] = {1,0,0};
 678   Double_t cosangleTorusSeg = (localVect2[0]*axisX[0]+
 679                                localVect2[1]*axisX[1]+
 680                                localVect2[2]*axisX[2])/normVect2;
 681   Double_t angleTorusSeg = TMath::ACos(cosangleTorusSeg)*TMath::RadToDeg();
 682
 683   TGeoRotation rotTorus("",angleRot1*TMath::RadToDeg(),
 684                         angleRotDiag*TMath::RadToDeg(),
 685                         180-angleTorusSeg+rotation);
 686   rotTorus.MultiplyBy(&rotTorusToZ, kTRUE);
 687   rotTorus.MasterToLocalVect(vect2, localVect2);
 688   if (localVect2[1]<0) {
 689     localVect2[0] = -localVect2[0];
 690     localVect2[1] = -localVect2[1];
 691     localVect2[2] = -localVect2[2];
 692   };
 693   normVect2 = TMath::Sqrt(localVect2[0]*localVect2[0]+
 694                           localVect2[1]*localVect2[1]+
 695                           localVect2[2]*localVect2[2]);
 696   Double_t axisY[3] = {0,1,0};
 697   Double_t cosPhi = (localVect2[0]*axisY[0]+localVect2[1]*axisY[1]+
 698                      localVect2[2]*axisY[2])/normVect2;
 699   Double_t torusPhi1 = TMath::ACos(cosPhi);
 700   Double_t torusR = (length/2)/TMath::Sin(torusPhi1);
 701   torusPhi1 = torusPhi1*TMath::RadToDeg();
 702   Double_t perpLength = TMath::Sqrt(torusR*torusR-length*length/4);
 703   Double_t localTransT[3] = {-perpLength,0,0};
 704   Double_t globalTransT[3];
 705   rotTorus.LocalToMasterVect(localTransT, globalTransT);
 706   TGeoTranslation transTorus("",cx+globalTransT[0],cy+globalTransT[1],
 707                              cz+globalTransT[2]);
 708
 709   TGeoCombiTrans  *combiTorus = new TGeoCombiTrans(transTorus, rotTorus);
 710
 711   //=================================================
 712   // Create the segment and add it to the mother volume
 713   TGeoVolume *vCableSegT = CreateTorus(torusPhi1, torusR, p2);
 714   p2Vol->AddNode(vCableSegT, p2, combiTorus);
 715
 716   if (fDebug) {
 717     printf("---\n  Cable segment points : ");
 718     printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
 719     printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
 720   };
 721
 722   if (ct) *ct = combiTorus;
 723   return vCableSegT;
 724 }
 725
 726 //________________________________________________________________________
 727 TGeoVolume *AliITSv11GeomCableRound::CreateSegment( Double_t *coord1,
 728                                                       Double_t *coord2,
 729                                                       Double_t *localVect1,
 730                                                       Double_t *localVect2, Int_t p)
 731 {
 732   // Create a cylindrical segment and its layers. The tube section is cutted by
 733   // two planes, defined by the normal vectors localVect1 and localVect2
 734
 735   //=================================================
 736   // Calculate segment "deformation"
 737   Double_t dx = coord2[0]-coord1[0];
 738   Double_t dy = coord2[1]-coord1[1];
 739   Double_t dz = coord2[2]-coord1[2];
 740   Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
 741
 742   // normal vectors have to point outside the TGeoCtub :
 743   if (-localVect1[2]<0) {
 744     localVect1[0] = -localVect1[0];
 745     localVect1[1] = -localVect1[1];
 746     localVect1[2] = -localVect1[2];
 747   };
 748   if (localVect2[2]<0) {
 749     localVect2[0] = -localVect2[0];
 750     localVect2[1] = -localVect2[1];
 751     localVect2[2] = -localVect2[2];
 752   };
 753   //=================================================
 754   // Create the segment
 755   TGeoCtub *cableSeg = new TGeoCtub(0, fRadius, length/2, fPhiMin, fPhiMax,
 756                                     localVect1[0],localVect1[1],localVect1[2],
 757                                     localVect2[0],localVect2[1],localVect2[2]);
 758
 759   TGeoMedium *skinMedia = fLayMedia[fNlayer-1];
 760   char name[100];
 761   sprintf(name, "%s_%i",GetName(), p);
 762   TGeoVolume *vCableSeg = new TGeoVolume(name, cableSeg, skinMedia);
 763   vCableSeg->SetLineColor(fLayColor[fNlayer-1]);
 764
 765   // add all cable layers
 766   Double_t layThickness[100+1];                        // 100 layers max !!!
 767   layThickness[0] = 0;
 768   for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
 769
 770     layThickness[iLay+1] = fLayThickness[iLay]+layThickness[iLay];
 771     TGeoCtub *lay = new TGeoCtub(layThickness[iLay], layThickness[iLay+1],
 772                                  length/2, fPhiMin, fPhiMax,
 773                                  localVect1[0],localVect1[1],localVect1[2],
 774                                  localVect2[0],localVect2[1],localVect2[2]);
 775
 776     TGeoVolume *vLay = new TGeoVolume("vCableSegLay", lay, fLayMedia[iLay]);
 777     vLay->SetLineColor(fLayColor[iLay]);
 778     vCableSeg->AddNode(vLay, iLay+1, 0);
 779   };
 780
 781   //vCableSeg->SetVisibility(kFALSE);
 782   return vCableSeg;
 783 }
 784
 785
 786 //________________________________________________________________________
 787 TGeoVolume *AliITSv11GeomCableRound::CreateTubeSegment( Double_t *coord1,
 788                                                       Double_t *coord2, Int_t p)
 789 {
 790   // Create a cylindrical segment and its layers
 791
 792   //=================================================
 793   // Calculate segment "deformation"
 794   Double_t dx = coord2[0]-coord1[0];
 795   Double_t dy = coord2[1]-coord1[1];
 796   Double_t dz = coord2[2]-coord1[2];
 797   Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
 798
 799  //=================================================
 800   // Create the segment
 801
 802   TGeoTubeSeg *cableSeg = new TGeoTubeSeg(0, fRadius, length/2, fPhiMin, fPhiMax);
 803
 804   TGeoMedium *skinMedia = fLayMedia[fNlayer-1];
 805   char name[100];
 806   sprintf(name, "%s_%i",GetName(), p);
 807   TGeoVolume *vCableSeg = new TGeoVolume(name, cableSeg, skinMedia);
 808   vCableSeg->SetLineColor(fLayColor[fNlayer-1]);
 809
 810   // add all cable layers
 811   Double_t layThickness[100+1];                        // 100 layers max !!!
 812   layThickness[0] = 0;
 813   for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
 814
 815     layThickness[iLay+1] = fLayThickness[iLay]+layThickness[iLay];
 816     TGeoTubeSeg*lay = new TGeoTubeSeg(layThickness[iLay], layThickness[iLay+1],
 817                                       length/2, fPhiMin, fPhiMax);
 818     TGeoVolume *vLay = new TGeoVolume("vCableSegLay", lay, fLayMedia[iLay]);
 819     vLay->SetLineColor(fLayColor[iLay]);
 820     vCableSeg->AddNode(vLay, iLay+1, 0);
 821   };
 822
 823   //vCableSeg->SetVisibility(kFALSE);
 824   return vCableSeg;
 825 }
 826
 827
 828 //________________________________________________________________________
 829 TGeoVolume *AliITSv11GeomCableRound::CreateTorus( Double_t &phi,
 830                                                   Double_t &r, Int_t p)
 831 {
 832   // Create one torus segment and its layers
 833
 834   Double_t torusR = r;
 835 //   Double_t torusPhi1 = phi;
 836 //   Double_t torusDPhi = -2*torusPhi1;  // bug in root ...
 837   Double_t torusPhi1 = 360-phi;
 838   Double_t torusDPhi = 2*phi;
 839
 840   //  // Create the segment, it will also work as the last layer
 841   TGeoTorus *cableSeg = new TGeoTorus(torusR, 0, fRadius, torusPhi1, torusDPhi);
 842   TGeoMedium *skinMedia = fLayMedia[fNlayer-1];
 843   char name[100];
 844   sprintf(name, "%s_%i",GetName(),p);
 845   TGeoVolume *vCableSeg = new TGeoVolume(name, cableSeg, skinMedia);
 846   vCableSeg->SetLineColor(fLayColor[fNlayer-1]);
 847
 848   // add all cable layers but last
 849   Double_t layThickness[100+1];                        // 100 layers max !!!
 850   layThickness[0] = 0;
 851   for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
 852
 853     layThickness[iLay+1] = fLayThickness[iLay]+layThickness[iLay];
 854     TGeoTorus *lay = new TGeoTorus(torusR, layThickness[iLay],
 855                                    layThickness[iLay+1],
 856                                    torusPhi1,torusDPhi);
 857
 858     TGeoVolume *vLay = new TGeoVolume("vCableSegLay",lay,fLayMedia[iLay]);
 859     vLay->SetLineColor(fLayColor[iLay]);
 860
 861     vCableSeg->AddNode(vLay, iLay+1,0);
 862   };
 863
 864   //vCableSeg->SetVisibility(kFALSE);
 865   return vCableSeg;
 866 }
 867
 868 //________________________________________________________________________
 869 void AliITSv11GeomCableRound::SetNLayers(Int_t nLayers) {
 870   // Set the total number of layers
 871   if((nLayers>0) &&(nLayers<=fgkCableMaxLayer)) {
 872     fNlayer = nLayers;
 873     for (Int_t i = 0; i<fNlayer; i++) {
 874       fLayThickness[i] = 0;
 875       fLayMedia[i] = 0;
 876     };
 877   };
 878 }
 879
 880 //________________________________________________________________________
 881 Int_t AliITSv11GeomCableRound::SetLayer(Int_t nLayer, Double_t thick,
 882                                            TGeoMedium *medium, Int_t color) {
 883   // Set layer #nLayer
 884   if ((nLayer<0)||(nLayer>=fNlayer)) {
 885     printf("Set wrong layer number of the cable\n");
 886     return kFALSE;
 887   };
 888   if (nLayer>0)
 889     if (fLayThickness[nLayer-1]<=0) {
 890       printf("You must define cable layer %i first !",nLayer-1);
 891       return kFALSE;
 892     };
 893
 894   Double_t thickTot = 0;
 895   for (Int_t i=0; i<nLayer; i++) thickTot += fLayThickness[i];
 896   thickTot += thick;
 897   if (thickTot-1e-10>fRadius) {
 898     printf("Can't add this layer, cable thickness would be higher than total\n");
 899     return kFALSE;
 900   };
 901
 902   fLayThickness[nLayer] = thick;
 903   fLayMedia[nLayer] = medium;
 904   fLayColor[nLayer] = color;
 905
 906   return kTRUE;
 907 }