]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ITS/AliITSv11GeomCableFlat.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Possibility of selecting one subdetector when making recpoints starting from raw...
[u/mrichter/AliRoot.git] / ITS / AliITSv11GeomCableFlat.cxx
CommitLineData
b7943f00 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// Class for flat cables
19//
20// Ludovic Gaudichet gaudichet@to.infn.it
21//*************************************************************************
22
23
24
25// General Root includes
26//#include <Riostream.h>
27#include <TMath.h>
28#include <TVectorD.h>
29
30// Root Geometry includes
31#include <TGeoManager.h>
32#include <TGeoVolume.h>
33#include <TGeoArb8.h>
7a82301d 34#include <TGeoTube.h>
b7943f00 35#include <TGeoMatrix.h>
36#include <TGeoNode.h>
37
38#include "AliITSv11GeomCableFlat.h"
39
40
41ClassImp(AliITSv11GeomCableFlat)
42
43//________________________________________________________________________
33ddec7d 44AliITSv11GeomCableFlat::AliITSv11GeomCableFlat():
45 AliITSv11GeomCable(),
46 fWidth(0),
47 fThick(0),
48 fNlayer(0)
b7943f00 49{
50 // constructor
b7943f00 51 for (Int_t i=0; i<fgkCableMaxLayer ; i++) {
52 fLayThickness[i] = 0;
53 fTranslation[i] = 0;
54 fLayColor[i] = 0;
55 fLayMedia[i] = 0;
56 };
5d7a6c6d 57}
b7943f00 58
59//________________________________________________________________________
60AliITSv11GeomCableFlat::
61AliITSv11GeomCableFlat(const char* name, Double_t width, Double_t thick) :
33ddec7d 62 AliITSv11GeomCable(name),
63 fWidth(width),
64 fThick(thick),
65 fNlayer(0)
66 {
b7943f00 67 // standard constructor
b7943f00 68 for (Int_t i=0; i<fgkCableMaxLayer ; i++) {
69 fLayThickness[i] = 0;
70 fTranslation[i] = 0;
71 fLayColor[i] = 0;
72 fLayMedia[i] = 0;
73 };
5d7a6c6d 74}
b7943f00 75
76//________________________________________________________________________
77AliITSv11GeomCableFlat::AliITSv11GeomCableFlat(const AliITSv11GeomCableFlat &s) :
78 AliITSv11GeomCable(s),fWidth(s.fWidth),fThick(s.fThick),fNlayer(s.fNlayer)
79{
80 // Copy Constructor
81 for (Int_t i=0; i<s.fNlayer; i++) {
82 fLayThickness[i] = s.fLayThickness[i];
83 fTranslation[i] = s.fTranslation[i];
84 fLayMedia[i] = s.fLayMedia[i];
85 fLayColor[i] = s.fLayColor[i];
86 }
87}
88
89//________________________________________________________________________
90AliITSv11GeomCableFlat& AliITSv11GeomCableFlat::
91operator=(const AliITSv11GeomCableFlat &s) {
92 // Assignment operator
93 // Not fully inplemented yet !!!
94
95 if(&s == this) return *this;
96 *this = s;
97 fWidth = s.fWidth;
98 fThick = s.fThick;
99 fNlayer = s.fNlayer;
100 for (Int_t i=0; i<s.fNlayer; i++) {
101 fLayThickness[i] = s.fLayThickness[i];
102 fTranslation[i] = s.fTranslation[i];
103 fLayMedia[i] = s.fLayMedia[i];
104 fLayColor[i] = s.fLayColor[i];
105 };
106 return *this;
107}
108
109//________________________________________________________________________
110Int_t AliITSv11GeomCableFlat::GetPoint( Int_t iCheckPt, Double_t *coord)
111 const {
112 // Get the correct point #iCheckPt
113 TVectorD *coordVector =(TVectorD *)fPointArray.At(2*iCheckPt);
114 if (coordVector) {
5d7a6c6d 115#if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
116 CopyFrom(coord, coordVector->GetElements());
117#else
b7943f00 118 CopyFrom(coord, coordVector->GetMatrixArray());
5d7a6c6d 119#endif
b7943f00 120 return kTRUE;
121 } else {
122 return kFALSE;
123 };
5d7a6c6d 124}
b7943f00 125
b7943f00 126//________________________________________________________________________
127Int_t AliITSv11GeomCableFlat::GetVect( Int_t iCheckPt, Double_t *coord)
128 const {
129 // Get the correct vect corresponding to point #iCheckPt
130
131 TVectorD *coordVector =(TVectorD *)fPointArray.At(2*iCheckPt+1);
132 if (coordVector) {
5d7a6c6d 133#if ROOT_VERSION_CODE < ROOT_VERSION(4,0,0)
134 CopyFrom(coord, coordVector->GetElements());
135#else
b7943f00 136 CopyFrom(coord, coordVector->GetMatrixArray());
5d7a6c6d 137#endif
b7943f00 138 return kTRUE;
139 } else {
140 return kFALSE;
141 };
5d7a6c6d 142}
b7943f00 143
b7943f00 144//________________________________________________________________________
145void AliITSv11GeomCableFlat::AddCheckPoint( TGeoVolume *vol, Int_t iCheckPt,
146 Double_t *coord, Double_t *orthVect)
147{
148 //
149 // Add a check point. In the fPointArray, the point is at i and its vector
150 // is at i+1.
151 //
152
153// if (iCheckPt>=fVolumeArray.GetEntriesFast()) {
154// fVolumeArray.AddLast(vol);
155// TVectorD *point = new TVectorD(3,coord);
156// TVectorD *vect = new TVectorD(3,orthVect);
157// fPointArray.AddLast(point);
158// fPointArray.AddLast(vect);
159
160// } else if ((iCheckPt >= 0)&&(iCheckPt < fVolumeArray.GetEntriesFast())) {
161// fVolumeArray.AddAt(vol, iCheckPt);
162// TVectorD *point = new TVectorD(3,coord);
163// TVectorD *vect = new TVectorD(3,orthVect);
164// fPointArray.AddAt(point, iCheckPt*2 );
165// fPointArray.AddAt(vect, iCheckPt*2+1);
166// };
167 fVolumeArray.AddAtAndExpand(vol, iCheckPt);
168 TVectorD *point = new TVectorD(3,coord);
169 TVectorD *vect = new TVectorD(3,orthVect);
170 fPointArray.AddAtAndExpand(point, iCheckPt*2 );
171 fPointArray.AddAtAndExpand(vect, iCheckPt*2+1);
5d7a6c6d 172}
b7943f00 173
174//________________________________________________________________________
175void AliITSv11GeomCableFlat::PrintCheckPoints() const {
176 // print all check points of the cable
177 printf(" ---\n Printing all check points of the flat cable\n");
178 for (Int_t i = 0; i<fVolumeArray.GetEntriesFast(); i++) {
179 Double_t coord[3];
180 if (GetPoint( i, coord))
181 printf(" ( %.2f, %.2f, %.2f )\n", coord[0], coord[1], coord[2]);
182 };
5d7a6c6d 183}
b7943f00 184
185//________________________________________________________________________
7a82301d 186TGeoVolume* AliITSv11GeomCableFlat::CreateAndInsertCableSegment(Int_t p2,
108bd0fe 187 Double_t rotation,
188 TGeoCombiTrans** ct)
b7943f00 189{
190// Creates a cable segment between points p1 and p2.
191// Rotation is the eventual rotation of the flat cable
192// along its length axis
193//
194// The segment volume is created inside the volume containing point2
195// Therefore this segment should be defined in this volume only.
196// I mean here that, if the previous point is in another volume,
197// it should be just at the border between the 2 volumes. Also the
198// orientation vector of the previous point should be orthogonal to
199// the surface between the 2 volumes.
200
201 TGeoNode *mainNode;
202 if (fInitialNode==0) {
203 TObjArray *nodes = gGeoManager->GetListOfNodes();
7a82301d 204 if (nodes->GetEntriesFast()==0) return 0;
b7943f00 205 mainNode = (TGeoNode *) nodes->UncheckedAt(0);
206 } else {
207 mainNode = fInitialNode;
208 };
209
210 Int_t p1 = p2 - 1;
211 TGeoVolume *p2Vol = GetVolume(p2);
212 TGeoVolume *p1Vol = GetVolume(p1);
213
214 ResetCheckDaughter();
215 fCurrentVol = p1Vol;
216 if (! CheckDaughter(mainNode)) {
217 printf("Error::volume containing point is not visible in node tree!\n");
7a82301d 218 return 0;
b7943f00 219 };
220
221 Double_t coord1[3], coord2[3], vect1[3], vect2[3];
222 //=================================================
223 // Get p1 position in the systeme of p2
224 if (p1Vol!=p2Vol) {
225
226 Int_t p1nodeInd[fgkCableMaxNodeLevel];
227 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
228 Int_t p1volLevel = 0;
229 while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
230 p1volLevel--;
231
232 ResetCheckDaughter();
233 fCurrentVol = p2Vol;
234 if (! CheckDaughter(mainNode)) {
235 printf("Error::volume containing point is not visible in node tree!\n");
7a82301d 236 return 0;
b7943f00 237 };
238 Int_t p2nodeInd[fgkCableMaxNodeLevel];
239 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
240 Int_t commonMotherLevel = 0;
241 while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
242 commonMotherLevel++;
243 commonMotherLevel--;
244 Int_t p2volLevel = 0;
245 while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
246 p2volLevel--;
247
248 // Get coord and vect of p1 in the common mother reference system
249 if (! GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1) )
7a82301d 250 return 0;
b7943f00 251 if (! GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1) )
7a82301d 252 return 0;
b7943f00 253
254 // Translate them in the reference system of the volume containing p2
255 TGeoNode *pathNode[fgkCableMaxNodeLevel];
256 pathNode[0] = mainNode;
257 for (Int_t i=0; i<=p2volLevel; i++) {
258 pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
259 };
260 Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
261 Double_t globalVect1[3] = {vect1[0], vect1[1], vect1[2]};
262
263 for (Int_t i = commonMotherLevel+1; i <= p2volLevel; i++) {
264 pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
265 pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
266 CopyFrom(globalCoord1, coord1);
267 CopyFrom(globalVect1, vect1);
268 };
269 } else {
7a82301d 270 if (! GetCheckPoint(p1, 0, 0, coord1) ) return 0;
271 if (! GetCheckVect(p1, 0, 0, vect1) ) return 0;
b7943f00 272 };
273
274 //=================================================
275 // Get p2 position in the systeme of p2
7a82301d 276 if (! GetCheckPoint(p2, 0, 0, coord2) ) return 0;
277 if (! GetCheckVect(p2, 0, 0, vect2) ) return 0;
b7943f00 278
279 Double_t cx = (coord1[0]+coord2[0])/2;
280 Double_t cy = (coord1[1]+coord2[1])/2;
281 Double_t cz = (coord1[2]+coord2[2])/2;
282 Double_t dx = coord2[0]-coord1[0];
283 Double_t dy = coord2[1]-coord1[1];
284 Double_t dz = coord2[2]-coord1[2];
285
286 //=================================================
287 // Positionning of the segment between the 2 points
288 if (TMath::Abs(dy)<1e-231) dy = 1e-231;
289 if (TMath::Abs(dz)<1e-231) dz = 1e-231;
290 //Double_t angleRot1 = -TMath::ATan(dx/dy);
291 //Double_t planDiagL = -TMath::Sqrt(dy*dy+dx*dx);
292 //if (dy<0) planDiagL = -planDiagL;
293 //Double_t angleRotDiag = TMath::ATan(planDiagL/dz);
294
295 Double_t angleRot1 = -TMath::ATan2(dx,dy);
296 Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
297 Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
298 //--- (Calculate rotation of segment on the Z axis)
299 //-- Here I'm trying to calculate the rotation to be applied in
300 //-- order to match as closer as possible this segment and the
301 //-- previous one.
302 //-- It seems that some times it doesn't work ...
b7943f00 303 TGeoRotation rotTemp("",angleRot1*TMath::RadToDeg(),
304 angleRotDiag*TMath::RadToDeg(), rotation);
305 Double_t localX[3] = {0,1,0};
306 Double_t globalX[3];
307 rotTemp.LocalToMasterVect(localX, globalX);
308 CopyFrom(localX, globalX);
309 GetCheckVect(localX, p2Vol, 0, fgkCableMaxNodeLevel+1, globalX);
310 Double_t orthVect[3];
311 GetCheckVect(vect1, p2Vol, 0, fgkCableMaxNodeLevel+1, orthVect);
73dfc864 312// Double_t angleRotZ = 0;
313// if (p2>1) {
314// Double_t orthVectNorm2 = ScalProd(orthVect,orthVect);
315// Double_t alpha1 = ScalProd(fPreviousX,orthVect)/orthVectNorm2;
316// Double_t alpha2 = ScalProd(globalX,orthVect)/orthVectNorm2;
317// Double_t globalX1p[3], globalX2p[3];
318// globalX1p[0] = fPreviousX[0] - alpha1*orthVect[0];
319// globalX1p[1] = fPreviousX[1] - alpha1*orthVect[1];
320// globalX1p[2] = fPreviousX[2] - alpha1*orthVect[2];
321// globalX2p[0] = globalX[0] - alpha2*orthVect[0];
322// globalX2p[1] = globalX[1] - alpha2*orthVect[1];
323// globalX2p[2] = globalX[2] - alpha2*orthVect[2];
324// //-- now I'm searching the 3th vect which makes an orthogonal base
325// //-- with orthVect and globalX1p ...
326// Double_t nulVect[3] = {0,0,0};
327// Double_t axis3[3];
328// TMath::Normal2Plane(nulVect, orthVect, globalX1p, axis3);
329// Double_t globalX1pNorm2 = ScalProd(globalX1p, globalX1p);
330// Double_t beta = ScalProd(globalX2p, globalX1p)/globalX1pNorm2;
331// Double_t gamma = ScalProd(globalX2p, axis3);
332// angleRotZ = (TMath::ATan2(1,0) - TMath::ATan2(beta, gamma))
333// *TMath::RadToDeg();
334// };
b7943f00 335 // cout << "!!!!!!!!!!!!!!!!!!! angle = " <<angleRotZ << endl;
336 CopyFrom(fPreviousX, globalX);
337 //---
338 Double_t localVect1[3], localVect2[3];
339 TGeoRotation rot("",angleRot1*TMath::RadToDeg(),
340 angleRotDiag*TMath::RadToDeg(),
341 rotation);
342// rotation-angleRotZ);
343// since angleRotZ doesn't always work, I won't use it ...
344
345 rot.MasterToLocalVect(vect1, localVect1);
346 rot.MasterToLocalVect(vect2, localVect2);
347
348 //=================================================
349 // Create the segment and add it to the mother volume
350 TGeoVolume *vCableSegB = CreateSegment(coord1, coord2,
351 localVect1, localVect2);
b7943f00 352 TGeoRotation rotArbSeg("", 0, 90, 0);
353 rotArbSeg.MultiplyBy(&rot, kFALSE);
354 TGeoTranslation trans("",cx, cy, cz);
355 TGeoCombiTrans *combiB = new TGeoCombiTrans(trans, rotArbSeg);
356 p2Vol->AddNode(vCableSegB, p2, combiB);
357 //=================================================;
358
359 if (fDebug) {
360 printf("---\n Cable segment points : ");
361 printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
362 printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
363 };
364
365// #include <TGeoSphere.h>
108bd0fe 366// TGeoMedium *airSDD = gGeoManager->GetMedium("ITS_AIR$");
b7943f00 367// TGeoSphere *sphere = new TGeoSphere(0, 0.05);
368// TGeoVolume *vSphere = new TGeoVolume("", sphere, airSDD);
369// TGeoTranslation *trC = new TGeoTranslation("", cx, cy, cz);
370// TGeoTranslation *tr1 = new TGeoTranslation("",coord1[0],
371// coord1[1],coord1[2]);
372// TGeoTranslation *tr2 = new TGeoTranslation("",coord2[0],
373// coord2[1],coord2[2]);
374// p2Vol->AddNode(vSphere, p2*3-2, trC);
375// p2Vol->AddNode(vSphere, p2*3-1, tr1);
376// p2Vol->AddNode(vSphere, p2*3 , tr2);
108bd0fe 377 if (ct) *ct = combiB;
7a82301d 378 return vCableSegB;
379}
380
381//________________________________________________________________________
382TGeoVolume* AliITSv11GeomCableFlat::CreateAndInsertBoxCableSegment(Int_t p2,
108bd0fe 383 Double_t rotation,
384 TGeoCombiTrans** ct)
7a82301d 385{
386 // This function is to be use only when the segment has the shape
387 // of a simple box, i.e. the normal vector to its end is perpendicular
388 // to the segment own axis
389// Creates a cable segment between points p1 and p2.
390// Rotation is the eventual rotation of the flat cable
391// along its length axis
392//
393// The segment volume is created inside the volume containing point2
394// Therefore this segment should be defined in this volume only.
395// I mean here that, if the previous point is in another volume,
396// it should be just at the border between the 2 volumes. Also the
397// orientation vector of the previous point should be orthogonal to
398// the surface between the 2 volumes.
399
400 TGeoNode *mainNode;
401 if (fInitialNode==0) {
402 TObjArray *nodes = gGeoManager->GetListOfNodes();
403 if (nodes->GetEntriesFast()==0) return 0;
404 mainNode = (TGeoNode *) nodes->UncheckedAt(0);
405 } else {
406 mainNode = fInitialNode;
407 };
408
409 Int_t p1 = p2 - 1;
410 TGeoVolume *p2Vol = GetVolume(p2);
411 TGeoVolume *p1Vol = GetVolume(p1);
412
413 ResetCheckDaughter();
414 fCurrentVol = p1Vol;
415 if (! CheckDaughter(mainNode)) {
416 printf("Error::volume containing point is not visible in node tree!\n");
417 return 0;
418 };
419
420 Double_t coord1[3], coord2[3], vect1[3], vect2[3];
421 //=================================================
422 // Get p1 position in the systeme of p2
423 if (p1Vol!=p2Vol) {
424
425 Int_t p1nodeInd[fgkCableMaxNodeLevel];
426 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
427 Int_t p1volLevel = 0;
428 while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
429 p1volLevel--;
430
431 ResetCheckDaughter();
432 fCurrentVol = p2Vol;
433 if (! CheckDaughter(mainNode)) {
434 printf("Error::volume containing point is not visible in node tree!\n");
435 return 0;
436 };
437 Int_t p2nodeInd[fgkCableMaxNodeLevel];
438 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
439 Int_t commonMotherLevel = 0;
440 while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
441 commonMotherLevel++;
442 commonMotherLevel--;
443 Int_t p2volLevel = 0;
444 while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
445 p2volLevel--;
446
447 // Get coord and vect of p1 in the common mother reference system
448 if (! GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1) )
449 return 0;
450 if (! GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1) )
451 return 0;
452
453 // Translate them in the reference system of the volume containing p2
454 TGeoNode *pathNode[fgkCableMaxNodeLevel];
455 pathNode[0] = mainNode;
456 for (Int_t i=0; i<=p2volLevel; i++) {
457 pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
458 };
459 Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
460 Double_t globalVect1[3] = {vect1[0], vect1[1], vect1[2]};
461
462 for (Int_t i = commonMotherLevel+1; i <= p2volLevel; i++) {
463 pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
464 pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
465 CopyFrom(globalCoord1, coord1);
466 CopyFrom(globalVect1, vect1);
467 };
468 } else {
469 if (! GetCheckPoint(p1, 0, 0, coord1) ) return 0;
470 if (! GetCheckVect(p1, 0, 0, vect1) ) return 0;
471 };
472
473 //=================================================
474 // Get p2 position in the systeme of p2
475 if (! GetCheckPoint(p2, 0, 0, coord2) ) return 0;
476 if (! GetCheckVect(p2, 0, 0, vect2) ) return 0;
477
478 Double_t cx = (coord1[0]+coord2[0])/2;
479 Double_t cy = (coord1[1]+coord2[1])/2;
480 Double_t cz = (coord1[2]+coord2[2])/2;
481 Double_t dx = coord2[0]-coord1[0];
482 Double_t dy = coord2[1]-coord1[1];
483 Double_t dz = coord2[2]-coord1[2];
484
485 //=================================================
486 // Positionning of the segment between the 2 points
487 if (TMath::Abs(dy)<1e-231) dy = 1e-231;
488 if (TMath::Abs(dz)<1e-231) dz = 1e-231;
489 //Double_t angleRot1 = -TMath::ATan(dx/dy);
490 //Double_t planDiagL = -TMath::Sqrt(dy*dy+dx*dx);
491 //if (dy<0) planDiagL = -planDiagL;
492 //Double_t angleRotDiag = TMath::ATan(planDiagL/dz);
493
494 Double_t angleRot1 = -TMath::ATan2(dx,dy);
495 Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
496 Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
497 //--- (Calculate rotation of segment on the Z axis)
498 //-- Here I'm trying to calculate the rotation to be applied in
499 //-- order to match as closer as possible this segment and the
500 //-- previous one.
501 //-- It seems that some times it doesn't work ...
7a82301d 502 TGeoRotation rotTemp("",angleRot1*TMath::RadToDeg(),
503 angleRotDiag*TMath::RadToDeg(), rotation);
504 Double_t localX[3] = {0,1,0};
505 Double_t globalX[3];
506 rotTemp.LocalToMasterVect(localX, globalX);
507 CopyFrom(localX, globalX);
508 GetCheckVect(localX, p2Vol, 0, fgkCableMaxNodeLevel+1, globalX);
509 Double_t orthVect[3];
510 GetCheckVect(vect1, p2Vol, 0, fgkCableMaxNodeLevel+1, orthVect);
73dfc864 511// Double_t angleRotZ = 0;
512// if (p2>1) {
513// Double_t orthVectNorm2 = ScalProd(orthVect,orthVect);
514// Double_t alpha1 = ScalProd(fPreviousX,orthVect)/orthVectNorm2;
515// Double_t alpha2 = ScalProd(globalX,orthVect)/orthVectNorm2;
516// Double_t globalX1p[3], globalX2p[3];
517// globalX1p[0] = fPreviousX[0] - alpha1*orthVect[0];
518// globalX1p[1] = fPreviousX[1] - alpha1*orthVect[1];
519// globalX1p[2] = fPreviousX[2] - alpha1*orthVect[2];
520// globalX2p[0] = globalX[0] - alpha2*orthVect[0];
521// globalX2p[1] = globalX[1] - alpha2*orthVect[1];
522// globalX2p[2] = globalX[2] - alpha2*orthVect[2];
523// //-- now I'm searching the 3th vect which makes an orthogonal base
524// //-- with orthVect and globalX1p ...
525// Double_t nulVect[3] = {0,0,0};
526// Double_t axis3[3];
527// TMath::Normal2Plane(nulVect, orthVect, globalX1p, axis3);
528// Double_t globalX1pNorm2 = ScalProd(globalX1p, globalX1p);
529// Double_t beta = ScalProd(globalX2p, globalX1p)/globalX1pNorm2;
530// Double_t gamma = ScalProd(globalX2p, axis3);
531// angleRotZ = (TMath::ATan2(1,0) - TMath::ATan2(beta, gamma))
532// *TMath::RadToDeg();
533// };
7a82301d 534 CopyFrom(fPreviousX, globalX);
535 //---
536 Double_t localVect1[3], localVect2[3];
537 TGeoRotation rot("",angleRot1*TMath::RadToDeg(),
538 angleRotDiag*TMath::RadToDeg(),
539 rotation);
540// rotation-angleRotZ);
541// since angleRotZ doesn't always work, I won't use it ...
542
543 rot.MasterToLocalVect(vect1, localVect1);
544 rot.MasterToLocalVect(vect2, localVect2);
545
546 //=================================================
547 // Create the segment and add it to the mother volume
548 TGeoVolume *vCableSegB = CreateBoxSegment(coord1, coord2);
549
550 TGeoRotation rotArbSeg("", 0, 90, 0);
551 rotArbSeg.MultiplyBy(&rot, kFALSE);
552 TGeoTranslation trans("",cx, cy, cz);
553 TGeoCombiTrans *combiB = new TGeoCombiTrans(trans, rotArbSeg);
554 p2Vol->AddNode(vCableSegB, p2, combiB);
555 //=================================================;
556
557 if (fDebug) {
558 printf("---\n Cable segment points : ");
559 printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
560 printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
561 };
562
108bd0fe 563 if (ct) *ct = combiB;
7a82301d 564 return vCableSegB;
565}
566
567//________________________________________________________________________
568TGeoVolume* AliITSv11GeomCableFlat::CreateAndInsertCableCylSegment(Int_t p2,
108bd0fe 569 Double_t rotation,
570 TGeoCombiTrans** ct)
7a82301d 571{
572 // Create a flat cable segment with a curvature between points p1 and p2.
573 // The radius and position of the curve is defined by the
574 // perpendicular vector of point p2 (the orientation of this vector
575 // and the position of the 2 check points are enough to completely
576 // define the curve)
577 // Rotation is the eventual rotation of the flat cable
578 // along its length axis
579 //
580
581 TGeoNode *mainNode;
582 if (fInitialNode==0) {
583 TObjArray *nodes = gGeoManager->GetListOfNodes();
584 if (nodes->GetEntriesFast()==0) return 0;
585 mainNode = (TGeoNode *) nodes->UncheckedAt(0);
586 } else {
587 mainNode = fInitialNode;
588 };
589
590 Int_t p1 = p2 - 1;
591 TGeoVolume *p1Vol = GetVolume(p1);
592 TGeoVolume *p2Vol = GetVolume(p2);
593
594 ResetCheckDaughter();
595 fCurrentVol = p1Vol;
596 if (! CheckDaughter(mainNode)) {
597 printf("Error::volume containing point is not visible in node tree!\n");
598 return 0;
599 };
600
601 Double_t coord1[3], coord2[3], vect1[3], vect2[3];
602 //=================================================
603 // Get p1 position in the systeme of p2
604 if (p1Vol!=p2Vol) {
605
606 Int_t p1nodeInd[fgkCableMaxNodeLevel];
607 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p1nodeInd[i]=fNodeInd[i];
608 Int_t p1volLevel = 0;
609 while (p1nodeInd[p1volLevel]!=-1) p1volLevel++;
610 p1volLevel--;
611
612 ResetCheckDaughter();
613 fCurrentVol = p2Vol;
614 if (! CheckDaughter(mainNode)) {
615 printf("Error::volume containing point is not visible in node tree!\n");
616 return 0;
617 };
618 Int_t p2nodeInd[fgkCableMaxNodeLevel];
619 for (Int_t i=0; i<fgkCableMaxNodeLevel; i++) p2nodeInd[i]=fNodeInd[i];
620 Int_t commonMotherLevel = 0;
621 while (p1nodeInd[commonMotherLevel]==fNodeInd[commonMotherLevel])
622 commonMotherLevel++;
623 commonMotherLevel--;
624 Int_t p2volLevel = 0;
625 while (fNodeInd[p2volLevel]!=-1) p2volLevel++;
626 p2volLevel--;
627
628 // Get coord and vect of p1 in the common mother reference system
629 GetCheckPoint(p1, 0, p1volLevel-commonMotherLevel, coord1);
630 GetCheckVect( p1, 0, p1volLevel-commonMotherLevel, vect1);
631 // Translate them in the reference system of the volume containing p2
632 TGeoNode *pathNode[fgkCableMaxNodeLevel];
633 pathNode[0] = mainNode;
634 for (Int_t i=0; i<=p2volLevel; i++) {
635 pathNode[i+1] = pathNode[i]->GetDaughter(p2nodeInd[i]);
636 };
637 Double_t globalCoord1[3] = {coord1[0], coord1[1], coord1[2]};
638 Double_t globalVect1[3] = {vect1[0], vect1[1], vect1[2]};
639
640 for (Int_t i = commonMotherLevel+1; i<=p2volLevel; i++) {
641 pathNode[i+1]->GetMatrix()->MasterToLocal(globalCoord1, coord1);
642 pathNode[i+1]->GetMatrix()->MasterToLocalVect(globalVect1, vect1);
643 CopyFrom(globalCoord1, coord1);
644 CopyFrom(globalVect1, vect1);
645 };
646 } else {
647 GetCheckPoint(p1, 0, 0, coord1);
648 GetCheckVect(p1, 0, 0, vect1);
649 };
650
651 //=================================================
652 // Get p2 position in the systeme of p2
653 GetCheckPoint(p2, 0, 0, coord2);
654 GetCheckVect(p2, 0, 0, vect2);
655
656 Double_t cx = (coord1[0]+coord2[0])/2;
657 Double_t cy = (coord1[1]+coord2[1])/2;
658 Double_t cz = (coord1[2]+coord2[2])/2;
659 Double_t dx = coord2[0]-coord1[0];
660 Double_t dy = coord2[1]-coord1[1];
661 Double_t dz = coord2[2]-coord1[2];
662 Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
663
664 //=================================================
665 // Positionning of the segment between the 2 points
666 if ((dy<1e-31)&&(dy>0)) dy = 1e-31;
667 if ((dz<1e-31)&&(dz>0)) dz = 1e-31;
668 if ((dy>-1e-31)&&(dy<0)) dy = -1e-31;
669 if ((dz>-1e-31)&&(dz<0)) dz = -1e-31;
670
671 Double_t angleRot1 = -TMath::ATan2(dx,dy);
672 Double_t planDiagL = TMath::Sqrt(dy*dy+dx*dx);
673 Double_t angleRotDiag = -TMath::ATan2(planDiagL,dz);
674
675 TGeoRotation rotTorusTemp("",angleRot1*TMath::RadToDeg(),
676 angleRotDiag*TMath::RadToDeg(),0);
677 TGeoRotation rotTorusToZ("",0,90,0);
678 rotTorusTemp.MultiplyBy(&rotTorusToZ, kTRUE);
679 Double_t localVect2[3];
680 rotTorusTemp.MasterToLocalVect(vect2, localVect2);
681 if (localVect2[1]<0) {
682 localVect2[0] = -localVect2[0];
683 localVect2[1] = -localVect2[1];
684 localVect2[2] = -localVect2[2];
685 };
686 Double_t normVect2 = TMath::Sqrt(localVect2[0]*localVect2[0]+
687 localVect2[1]*localVect2[1]+
688 localVect2[2]*localVect2[2]);
689 Double_t axisX[3] = {1,0,0};
690 Double_t cosangleTorusSeg = (localVect2[0]*axisX[0]+
691 localVect2[1]*axisX[1]+
692 localVect2[2]*axisX[2])/normVect2;
693 Double_t angleTorusSeg = TMath::ACos(cosangleTorusSeg)*TMath::RadToDeg();
694 TGeoRotation rotTorus("",angleRot1*TMath::RadToDeg(),
695 angleRotDiag*TMath::RadToDeg(),
696 45-angleTorusSeg+rotation);
697 //180-angleTorusSeg+rotation);
698 rotTorus.MultiplyBy(&rotTorusToZ, kTRUE);
699 rotTorus.MasterToLocalVect(vect2, localVect2);
700 if (localVect2[1]<0) {
701 localVect2[0] = -localVect2[0];
702 localVect2[1] = -localVect2[1];
703 localVect2[2] = -localVect2[2];
704 };
705 normVect2 = TMath::Sqrt(localVect2[0]*localVect2[0]+
706 localVect2[1]*localVect2[1]+
707 localVect2[2]*localVect2[2]);
708 Double_t axisY[3] = {0,1,0};
709 Double_t cosPhi = (localVect2[0]*axisY[0]+localVect2[1]*axisY[1]+
710 localVect2[2]*axisY[2])/normVect2;
711 Double_t torusPhi1 = TMath::ACos(cosPhi);
712 Double_t torusR = (length/2)/TMath::Sin(torusPhi1);
713 torusPhi1 = torusPhi1*TMath::RadToDeg();
714 Double_t perpLength = TMath::Sqrt(torusR*torusR-length*length/4);
715 Double_t localTransT[3] = {-perpLength,0,0};
716 Double_t globalTransT[3];
717 rotTorus.LocalToMasterVect(localTransT, globalTransT);
718 TGeoTranslation transTorus("",cx+globalTransT[0],cy+globalTransT[1],
719 cz+globalTransT[2]);
720
721 TGeoCombiTrans *combiTorus = new TGeoCombiTrans(transTorus, rotTorus);
722
723 //=================================================
724 // Create the segment and add it to the mother volume
725 TGeoVolume *vCableSegT = CreateCylSegment(torusPhi1, torusR);
726 p2Vol->AddNode(vCableSegT, p2, combiTorus);
727
728 if (fDebug) {
729 printf("---\n Cable segment points : ");
730 printf("%f, %f, %f\n",coord1[0], coord1[1], coord1[2]);
731 printf("%f, %f, %f\n",coord2[0], coord2[1], coord2[2]);
732 };
733
108bd0fe 734 if (ct) *ct = combiTorus;
7a82301d 735 return vCableSegT;
5d7a6c6d 736}
b7943f00 737
738//________________________________________________________________________
739TGeoVolume *AliITSv11GeomCableFlat::CreateSegment( Double_t *coord1,
7a82301d 740 Double_t *coord2,
741 Double_t *localVect1,
742 Double_t *localVect2 )
b7943f00 743{
73dfc864 744 // Create a segment with arbitrary vertices (general case)
b7943f00 745 //=================================================
746 // Calculate segment "deformation"
747 Double_t dx = coord2[0]-coord1[0];
748 Double_t dy = coord2[1]-coord1[1];
749 Double_t dz = coord2[2]-coord1[2];
750 Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
751
752 Double_t cosTheta1 = -1./TMath::Sqrt( 1 + localVect1[0]*localVect1[0]
753 /localVect1[2]/localVect1[2] );
754 Double_t cosTheta2 = 1./TMath::Sqrt( 1 + localVect2[0]*localVect2[0]
755 /localVect2[2]/localVect2[2] );
756 if (localVect1[2]<0) cosTheta1 = -cosTheta1;
757 if (localVect2[2]<0) cosTheta2 = -cosTheta2;
758
759 Double_t dL1 = 0.5*fWidth*TMath::Tan(TMath::ACos(cosTheta1));
760 Double_t dL2 = 0.5*fWidth*TMath::Tan(TMath::ACos(cosTheta2));
761 if (localVect1[0]<0) dL1 = - dL1;
762 if (localVect2[0]<0) dL2 = - dL2;
763 //---
764 Double_t cosPhi1 = -1./TMath::Sqrt( 1 + localVect1[1]*localVect1[1]
765 /localVect1[2]/localVect1[2] );
766 Double_t cosPhi2 = 1./TMath::Sqrt( 1 + localVect2[1]*localVect2[1]
767 /localVect2[2]/localVect2[2] );
768 if (localVect1[2]<0) cosPhi1 = -cosPhi1;
769 if (localVect2[2]<0) cosPhi2 = -cosPhi2;
770
771 Double_t tanACosCosPhi1 = TMath::Tan(TMath::ACos(cosPhi1));
772 Double_t tanACosCosPhi2 = TMath::Tan(TMath::ACos(cosPhi2));
773 if (localVect1[1]<0) tanACosCosPhi1 = -tanACosCosPhi1;
774 if (localVect2[1]<0) tanACosCosPhi2 = -tanACosCosPhi2;
775
108bd0fe 776 Double_t dl1 = 0.5*fThick*tanACosCosPhi1*0.99999999999999;
777 Double_t dl2 = 0.5*fThick*tanACosCosPhi2*0.99999999999999;
73dfc864 778 // 0.9999999999999 is for correcting problems in TGeo...
b7943f00 779 //=================================================
780 // Create the segment
781 TGeoArb8 *cableSeg = new TGeoArb8(fThick/2);
782 cableSeg->SetVertex( 0, -fWidth/2, -length/2 - dL1 + dl1);
531d6cdc 783 cableSeg->SetVertex( 1, -fWidth/2, length/2 + dL2 - dl2);
b7943f00 784 cableSeg->SetVertex( 2, fWidth/2, length/2 - dL2 - dl2);
531d6cdc 785 cableSeg->SetVertex( 3, fWidth/2, -length/2 + dL1 + dl1);
b7943f00 786 cableSeg->SetVertex( 4, -fWidth/2, -length/2 - dL1 - dl1);
531d6cdc 787 cableSeg->SetVertex( 5, -fWidth/2, length/2 + dL2 + dl2);
b7943f00 788 cableSeg->SetVertex( 6, fWidth/2, length/2 - dL2 + dl2);
531d6cdc 789 cableSeg->SetVertex( 7, fWidth/2, -length/2 + dL1 - dl1);
b7943f00 790
73dfc864 791 TGeoVolume *vCableSeg = new TGeoVolume(GetName(), cableSeg, fLayMedia[fNlayer-1]);
b7943f00 792
73dfc864 793 // add all cable layers but the last
794 for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
b7943f00 795
796 Double_t dl1Lay = 0.5*fLayThickness[iLay]*tanACosCosPhi1;
797 Double_t dl2Lay = 0.5*fLayThickness[iLay]*tanACosCosPhi2;
798
799 Double_t ztr = -fThick/2;
800 for (Int_t i=0;i<iLay; i++) ztr+= fLayThickness[i];
801 ztr+= fLayThickness[iLay]/2;
802
803 Double_t dl1LayS = ztr*tanACosCosPhi1;
804 Double_t dl2LayS = ztr*tanACosCosPhi2;
805
806 TGeoArb8 *lay = new TGeoArb8(fLayThickness[iLay]/2);
807 lay->SetVertex( 0, -fWidth/2, -length/2 - dL1 + dl1Lay - dl1LayS);
531d6cdc 808 lay->SetVertex( 1, -fWidth/2, length/2 + dL2 - dl2Lay + dl2LayS);
b7943f00 809 lay->SetVertex( 2, fWidth/2, length/2 - dL2 - dl2Lay + dl2LayS);
531d6cdc 810 lay->SetVertex( 3, fWidth/2, -length/2 + dL1 + dl1Lay - dl1LayS);
b7943f00 811 lay->SetVertex( 4, -fWidth/2, -length/2 - dL1 - dl1Lay - dl1LayS);
531d6cdc 812 lay->SetVertex( 5, -fWidth/2, length/2 + dL2 + dl2Lay + dl2LayS);
b7943f00 813 lay->SetVertex( 6, fWidth/2, length/2 - dL2 + dl2Lay + dl2LayS);
531d6cdc 814 lay->SetVertex( 7, fWidth/2, -length/2 + dL1 - dl1Lay - dl1LayS);
b7943f00 815 TGeoVolume *vLay = new TGeoVolume("vCableSegLay", lay, fLayMedia[iLay]);
816 vLay->SetLineColor(fLayColor[iLay]);
817
818 if (fTranslation[iLay]==0)
819 fTranslation[iLay] = new TGeoTranslation(0, 0, ztr);
820 vCableSeg->AddNode(vLay, iLay+1, fTranslation[iLay]);
821 };
822
73dfc864 823 //vCableSeg->SetVisibility(kFALSE);
7a82301d 824 return vCableSeg;
825}
826
7a82301d 827//________________________________________________________________________
828TGeoVolume *AliITSv11GeomCableFlat::CreateCylSegment(Double_t &phi,
829 Double_t &r)
830{
73dfc864 831 // Create a segment in shape of a cylinder, allows to represent
832 // a folded flat cable
7a82301d 833
834 Double_t phi1 = 360-phi;
835 Double_t phi2 = 360+phi;
836
837 Double_t rMin = r-fThick/2;
838 Double_t rMax = r+fThick/2;
839 //=================================================
840 // Create the segment
841
842 TGeoTubeSeg *cableSeg = new TGeoTubeSeg(rMin, rMax, fWidth/2,
843 phi1, phi2);
73dfc864 844 TGeoVolume *vCableSeg = new TGeoVolume(GetName(), cableSeg, fLayMedia[fNlayer-1]);
7a82301d 845
73dfc864 846 // add all cable layers but the last
847 for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
7a82301d 848
849 Double_t ztr = -fThick/2;
850 for (Int_t i=0;i<iLay; i++) ztr+= fLayThickness[i];
851
852 rMin = r + ztr;
853 rMax = r + ztr + fLayThickness[iLay];
854 TGeoTubeSeg *lay = new TGeoTubeSeg(rMin, rMax, fWidth/2,
855 phi1, phi2);
856
857 TGeoVolume *vLay = new TGeoVolume("vCableSegLay", lay, fLayMedia[iLay]);
858 vLay->SetLineColor(fLayColor[iLay]);
859
860 vCableSeg->AddNode(vLay, iLay+1, 0);
861 };
862
73dfc864 863 //vCableSeg->SetVisibility(kFALSE);
b7943f00 864 return vCableSeg;
5d7a6c6d 865}
b7943f00 866
7a82301d 867//________________________________________________________________________
868TGeoVolume *AliITSv11GeomCableFlat::CreateBoxSegment( Double_t *coord1,
869 Double_t *coord2)
870{
7a82301d 871 // Create a segment for the case it is a simple box
73dfc864 872 //=================================================
7a82301d 873 Double_t dx = coord2[0]-coord1[0];
874 Double_t dy = coord2[1]-coord1[1];
875 Double_t dz = coord2[2]-coord1[2];
876 Double_t length = TMath::Sqrt(dx*dx+dy*dy+dz*dz);
877
878 TGeoBBox *cableSeg = new TGeoBBox(fWidth/2, length/2, fThick/2);
73dfc864 879 TGeoVolume *vCableSeg = new TGeoVolume(GetName(), cableSeg, fLayMedia[fNlayer-1]);
880 // This volume is the cable container. It codes also the material for the
881 // last layer
7a82301d 882
73dfc864 883 // add all cable layers but the last one
884 for (Int_t iLay=0; iLay<fNlayer-1; iLay++) {
7a82301d 885
886 Double_t ztr = -fThick/2;
887 for (Int_t i=0;i<iLay; i++) ztr+= fLayThickness[i];
888 ztr+= fLayThickness[iLay]/2;
889
890 TGeoBBox *lay = new TGeoBBox(fWidth/2, length/2, fLayThickness[iLay]/2);
891
7a82301d 892 TGeoVolume *vLay = new TGeoVolume("vCableSegLay", lay, fLayMedia[iLay]);
893 vLay->SetLineColor(fLayColor[iLay]);
894
895 if (fTranslation[iLay]==0)
896 fTranslation[iLay] = new TGeoTranslation(0, 0, ztr);
897 vCableSeg->AddNode(vLay, iLay+1, fTranslation[iLay]);
898 };
899
73dfc864 900 //vCableSeg->SetVisibility(kFALSE);
7a82301d 901 return vCableSeg;
902}
903
b7943f00 904//________________________________________________________________________
905void AliITSv11GeomCableFlat::SetNLayers(Int_t nLayers) {
906 // Set the number of layers
907 if((nLayers>0) &&(nLayers<=fgkCableMaxLayer)) {
908
909 fNlayer = nLayers;
910 for (Int_t i=0; i<fgkCableMaxLayer ; i++) {
911 fLayThickness[i] = 0;
912 fTranslation[i] = 0;
913 fLayColor[i] = 0;
914 fLayMedia[i] = 0;
915 };
916 };
5d7a6c6d 917}
b7943f00 918
919//________________________________________________________________________
920Int_t AliITSv11GeomCableFlat::SetLayer(Int_t nLayer, Double_t thick,
921 TGeoMedium *medium, Int_t color) {
922 // Set the layer number nLayer
923 if ((nLayer<0)||(nLayer>=fNlayer)) {
924 printf("Set wrong layer number of the cable\n");
925 return kFALSE;
926 };
927 if (nLayer>0)
928 if (fLayThickness[nLayer-1]<=0) {
929 printf("AliITSv11GeomCableFlat::SetLayer():"
930 " You must define cable layer %i first !",nLayer-1);
931 return kFALSE;
932 };
933
934 Double_t thickTot = 0;
935 for (Int_t i=0; i<nLayer; i++) thickTot += fLayThickness[i];
936 thickTot += thick;
937 if (thickTot-1e-10>fThick) {
938 printf("Can't add this layer, cable thickness would be higher than total\n");
939 return kFALSE;
940 };
941
942 fLayThickness[nLayer] = thick;
943 fLayMedia[nLayer] = medium;
944 fLayColor[nLayer] = color;
945 fTranslation[nLayer] = 0;
946 return kTRUE;
5d7a6c6d 947}