Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoTube.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Andrei Gheata 24/10/01
3// TGeoTube::Contains() and DistFromInside/In() implemented by Mihaela Gheata
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13/** \class TGeoTube
14\ingroup Geometry_classes
15
16Cylindrical tube class. A tube has 3 parameters :
17 - Rmin - minimum radius
18 - Rmax - maximum radius
19 - dz - half length
20Begin_Macro(source)
21{
22 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
23 new TGeoManager("tube", "poza2");
24 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
25 TGeoMedium *med = new TGeoMedium("MED",1,mat);
26 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
27 gGeoManager->SetTopVolume(top);
28 TGeoVolume *vol = gGeoManager->MakeTube("TUBE",med, 20,30,40);
29 vol->SetLineWidth(2);
30 top->AddNode(vol,1);
31 gGeoManager->CloseGeometry();
32 gGeoManager->SetNsegments(80);
33 top->Draw();
34 TView *view = gPad->GetView();
35 view->ShowAxis();
36}
37End_Macro
38*/
39
40/** \class TGeoTubeSeg
41\ingroup Geometry_classes
42
43A phi segment of a tube. Has 5 parameters :
44 - the same 3 as a tube;
45 - first phi limit (in degrees)
46 - second phi limit
47The segment will be be placed from the first angle (first phi limit)
48to the second angle (second phi limit)
49Begin_Macro(source)
50{
51 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
52 new TGeoManager("tubeseg", "poza3");
53 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
54 TGeoMedium *med = new TGeoMedium("MED",1,mat);
55 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
56 gGeoManager->SetTopVolume(top);
57 TGeoVolume *vol = gGeoManager->MakeTubs("TUBESEG",med, 20,30,40,-30,270);
58 vol->SetLineWidth(2);
59 top->AddNode(vol,1);
60 gGeoManager->CloseGeometry();
61 gGeoManager->SetNsegments(80);
62 top->Draw();
63 TView *view = gPad->GetView();
64 view->ShowAxis();
65}
66
67End_Macro
68*/
69
70/** \class TGeoCtub
71\ingroup Geometry_classes
72
73A tube segment cut with 2 planes. Has 11 parameters :
74 - the same 5 as a tube segment;
75 - x, y, z components of the normal to the -dZ cut plane in
76 point (0, 0, -dZ);
77 - x, y, z components of the normal to the +dZ cut plane in
78 point (0, 0, dZ);
79
80Begin_Macro(source)
81{
82 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
83 new TGeoManager("ctub", "poza3");
84 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
85 TGeoMedium *med = new TGeoMedium("MED",1,mat);
86 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
87 gGeoManager->SetTopVolume(top);
88 Double_t theta = 160.*TMath::Pi()/180.;
89 Double_t phi = 30.*TMath::Pi()/180.;
90 Double_t nlow[3];
91 nlow[0] = TMath::Sin(theta)*TMath::Cos(phi);
92 nlow[1] = TMath::Sin(theta)*TMath::Sin(phi);
93 nlow[2] = TMath::Cos(theta);
94 theta = 20.*TMath::Pi()/180.;
95 phi = 60.*TMath::Pi()/180.;
96 Double_t nhi[3];
97 nhi[0] = TMath::Sin(theta)*TMath::Cos(phi);
98 nhi[1] = TMath::Sin(theta)*TMath::Sin(phi);
99 nhi[2] = TMath::Cos(theta);
100 TGeoVolume *vol = gGeoManager->MakeCtub("CTUB",med, 20,30,40,-30,250, nlow[0], nlow[1], nlow[2], nhi[0],nhi[1],nhi[2]);
101 vol->SetLineWidth(2);
102 top->AddNode(vol,1);
103 gGeoManager->CloseGeometry();
104 gGeoManager->SetNsegments(80);
105 top->Draw();
106 TView *view = gPad->GetView();
107 view->ShowAxis();
108}
109End_Macro
110*/
111
112#include <iostream>
113
114#include "TGeoManager.h"
115#include "TGeoVolume.h"
116#include "TVirtualGeoPainter.h"
117#include "TGeoTube.h"
118#include "TBuffer3D.h"
119#include "TBuffer3DTypes.h"
120#include "TMath.h"
121
123
124////////////////////////////////////////////////////////////////////////////////
125/// Default constructor
126
128{
130 fRmin = 0.0;
131 fRmax = 0.0;
132 fDz = 0.0;
133}
134
135
136////////////////////////////////////////////////////////////////////////////////
137/// Default constructor specifying minimum and maximum radius
138
140 :TGeoBBox(0, 0, 0)
141{
143 SetTubeDimensions(rmin, rmax, dz);
144 if ((fDz<0) || (fRmin<0) || (fRmax<0)) {
146// if (fRmax<=fRmin) SetShapeBit(kGeoInvalidShape);
147// printf("tube : dz=%f rmin=%f rmax=%f\n", dz, rmin, rmax);
148 }
149 ComputeBBox();
150}
151////////////////////////////////////////////////////////////////////////////////
152/// Default constructor specifying minimum and maximum radius
153
155 :TGeoBBox(name, 0, 0, 0)
156{
158 SetTubeDimensions(rmin, rmax, dz);
159 if ((fDz<0) || (fRmin<0) || (fRmax<0)) {
161// if (fRmax<=fRmin) SetShapeBit(kGeoInvalidShape);
162// printf("tube : dz=%f rmin=%f rmax=%f\n", dz, rmin, rmax);
163 }
164 ComputeBBox();
165}
166
167////////////////////////////////////////////////////////////////////////////////
168/// Default constructor specifying minimum and maximum radius
169/// - param[0] = Rmin
170/// - param[1] = Rmax
171/// - param[2] = dz
172
174 :TGeoBBox(0, 0, 0)
175{
177 SetDimensions(param);
178 if ((fDz<0) || (fRmin<0) || (fRmax<0)) SetShapeBit(kGeoRunTimeShape);
179 ComputeBBox();
180}
181
182////////////////////////////////////////////////////////////////////////////////
183/// destructor
184
186{
187}
188
189////////////////////////////////////////////////////////////////////////////////
190/// Computes capacity of the shape in [length^3]
191
193{
195}
196
197////////////////////////////////////////////////////////////////////////////////
198/// Computes capacity of the shape in [length^3]
199
201{
202 Double_t capacity = 2.*TMath::Pi()*(rmax*rmax-rmin*rmin)*dz;
203 return capacity;
204}
205
206////////////////////////////////////////////////////////////////////////////////
207/// compute bounding box of the tube
208
210{
211 fDX = fDY = fRmax;
212 fDZ = fDz;
213}
214
215////////////////////////////////////////////////////////////////////////////////
216/// Compute normal to closest surface from POINT.
217
218void TGeoTube::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
219{
220 Double_t saf[3];
221 Double_t rsq = point[0]*point[0]+point[1]*point[1];
222 Double_t r = TMath::Sqrt(rsq);
223 saf[0] = TMath::Abs(fDz-TMath::Abs(point[2]));
224 saf[1] = (fRmin>1E-10)?TMath::Abs(r-fRmin):TGeoShape::Big();
225 saf[2] = TMath::Abs(fRmax-r);
226 Int_t i = TMath::LocMin(3,saf);
227 if (i==0) {
228 norm[0] = norm[1] = 0.;
229 norm[2] = TMath::Sign(1.,dir[2]);
230 return;
231 }
232 norm[2] = 0;
233 Double_t phi = TMath::ATan2(point[1], point[0]);
234 norm[0] = TMath::Cos(phi);
235 norm[1] = TMath::Sin(phi);
236 if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
237 norm[0] = -norm[0];
238 norm[1] = -norm[1];
239 }
240}
241
242////////////////////////////////////////////////////////////////////////////////
243/// Compute normal to closest surface from POINT.
244
245void TGeoTube::ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm,
246 Double_t /*rmin*/, Double_t /*rmax*/, Double_t /*dz*/)
247{
248 norm[2] = 0;
249 Double_t phi = TMath::ATan2(point[1], point[0]);
250 norm[0] = TMath::Cos(phi);
251 norm[1] = TMath::Sin(phi);
252 if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
253 norm[0] = -norm[0];
254 norm[1] = -norm[1];
255 }
256}
257
258////////////////////////////////////////////////////////////////////////////////
259/// test if point is inside this tube
260
262{
263 if (TMath::Abs(point[2]) > fDz) return kFALSE;
264 Double_t r2 = point[0]*point[0]+point[1]*point[1];
265 if ((r2<fRmin*fRmin) || (r2>fRmax*fRmax)) return kFALSE;
266 return kTRUE;
267}
268
269////////////////////////////////////////////////////////////////////////////////
270/// compute closest distance from point px,py to each corner
271
273{
275 Int_t numPoints = 4*n;
276 if (!HasRmin()) numPoints = 2*(n+1);
277 return ShapeDistancetoPrimitive(numPoints, px, py);
278}
279
280////////////////////////////////////////////////////////////////////////////////
281/// Compute distance from inside point to surface of the tube (static)
282/// Boundary safe algorithm.
283/// compute distance to surface
284/// Do Z
285
287{
289 if (dir[2]) {
290 sz = (TMath::Sign(dz, dir[2])-point[2])/dir[2];
291 if (sz<=0) return 0.0;
292 }
293 // Do R
294 Double_t nsq=dir[0]*dir[0]+dir[1]*dir[1];
295 if (TMath::Abs(nsq)<TGeoShape::Tolerance()) return sz;
296 Double_t rsq=point[0]*point[0]+point[1]*point[1];
297 Double_t rdotn=point[0]*dir[0]+point[1]*dir[1];
298 Double_t b,d;
300 // inner cylinder
301 if (rmin>0) {
302 // Protection in case point is actually outside the tube
303 if (rsq <= rmin*rmin+TGeoShape::Tolerance()) {
304 if (rdotn<0) return 0.0;
305 } else {
306 if (rdotn<0) {
307 DistToTube(rsq,nsq,rdotn,rmin,b,d);
308 if (d>0) {
309 sr=-b-d;
310 if (sr>0) return TMath::Min(sz,sr);
311 }
312 }
313 }
314 }
315 // outer cylinder
316 if (rsq >= rmax*rmax-TGeoShape::Tolerance()) {
317 if (rdotn>=0) return 0.0;
318 }
319 DistToTube(rsq,nsq,rdotn,rmax,b,d);
320 if (d>0) {
321 sr=-b+d;
322 if (sr>0) return TMath::Min(sz,sr);
323 }
324 return 0.;
325}
326
327////////////////////////////////////////////////////////////////////////////////
328/// Compute distance from inside point to surface of the tube
329/// Boundary safe algorithm.
330
331Double_t TGeoTube::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
332{
333 if (iact<3 && safe) {
334 *safe = Safety(point, kTRUE);
335 if (iact==0) return TGeoShape::Big();
336 if ((iact==1) && (*safe>step)) return TGeoShape::Big();
337 }
338 // compute distance to surface
339 return DistFromInsideS(point, dir, fRmin, fRmax, fDz);
340}
341
342////////////////////////////////////////////////////////////////////////////////
343/// Static method to compute distance from outside point to a tube with given parameters
344/// Boundary safe algorithm.
345/// check Z planes
346
348{
349 Double_t xi,yi,zi;
350 Double_t rmaxsq = rmax*rmax;
351 Double_t rminsq = rmin*rmin;
352 zi = dz - TMath::Abs(point[2]);
354 Bool_t in = kFALSE;
355 Bool_t inz = (zi<0)?kFALSE:kTRUE;
356 if (!inz) {
357 if (point[2]*dir[2]>=0) return TGeoShape::Big();
358 s = -zi/TMath::Abs(dir[2]);
359 xi = point[0]+s*dir[0];
360 yi = point[1]+s*dir[1];
361 Double_t r2=xi*xi+yi*yi;
362 if ((rminsq<=r2) && (r2<=rmaxsq)) return s;
363 }
364
365 Double_t rsq = point[0]*point[0]+point[1]*point[1];
366 // check outer cyl. surface
367 Double_t nsq=dir[0]*dir[0]+dir[1]*dir[1];
368 Double_t rdotn=point[0]*dir[0]+point[1]*dir[1];
369 Double_t b,d;
370 Bool_t inrmax = kFALSE;
371 Bool_t inrmin = kFALSE;
372 if (rsq<=rmaxsq+TGeoShape::Tolerance()) inrmax = kTRUE;
373 if (rsq>=rminsq-TGeoShape::Tolerance()) inrmin = kTRUE;
374 in = inz & inrmin & inrmax;
375 // If inside, we are most likely on a boundary within machine precision.
376 if (in) {
377 Bool_t checkout = kFALSE;
378 Double_t r = TMath::Sqrt(rsq);
379 if (zi<rmax-r) {
380 if ((TGeoShape::IsSameWithinTolerance(rmin,0)) || (zi<r-rmin)) {
381 if (point[2]*dir[2]<0) return 0.0;
382 return TGeoShape::Big();
383 }
384 }
385 if ((rmaxsq-rsq) < (rsq-rminsq)) checkout = kTRUE;
386 if (checkout) {
387 if (rdotn>=0) return TGeoShape::Big();
388 return 0.0;
389 }
390 if (TGeoShape::IsSameWithinTolerance(rmin,0)) return 0.0;
391 if (rdotn>=0) return 0.0;
392 // Ray exiting rmin -> check (+) solution for inner tube
394 DistToTube(rsq, nsq, rdotn, rmin, b, d);
395 if (d>0) {
396 s=-b+d;
397 if (s>0) {
398 zi=point[2]+s*dir[2];
399 if (TMath::Abs(zi)<=dz) return s;
400 }
401 }
402 return TGeoShape::Big();
403 }
404 // Check outer cylinder (only r>rmax has to be considered)
406 if (!inrmax) {
407 DistToTube(rsq, nsq, rdotn, rmax, b, d);
408 if (d>0) {
409 s=-b-d;
410 if (s>0) {
411 zi=point[2]+s*dir[2];
412 if (TMath::Abs(zi)<=dz) return s;
413 }
414 }
415 }
416 // check inner cylinder
417 if (rmin>0) {
418 DistToTube(rsq, nsq, rdotn, rmin, b, d);
419 if (d>0) {
420 s=-b+d;
421 if (s>0) {
422 zi=point[2]+s*dir[2];
423 if (TMath::Abs(zi)<=dz) return s;
424 }
425 }
426 }
427 return TGeoShape::Big();
428}
429
430////////////////////////////////////////////////////////////////////////////////
431/// Compute distance from outside point to surface of the tube and safe distance
432/// Boundary safe algorithm.
433/// fist localize point w.r.t tube
434
435Double_t TGeoTube::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
436{
437 if (iact<3 && safe) {
438 *safe = Safety(point, kFALSE);
439 if (iact==0) return TGeoShape::Big();
440 if ((iact==1) && (step<=*safe)) return TGeoShape::Big();
441 }
442// Check if the bounding box is crossed within the requested distance
443 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
444 if (sdist>=step) return TGeoShape::Big();
445 // find distance to shape
446 return DistFromOutsideS(point, dir, fRmin, fRmax, fDz);
447}
448
449////////////////////////////////////////////////////////////////////////////////
450/// Static method computing the distance to a tube with given radius, starting from
451/// POINT along DIR director cosines. The distance is computed as :
452/// RSQ = point[0]*point[0]+point[1]*point[1]
453/// NSQ = dir[0]*dir[0]+dir[1]*dir[1] ---> should NOT be 0 !!!
454/// RDOTN = point[0]*dir[0]+point[1]*dir[1]
455/// The distance can be computed as :
456/// D = -B +/- DELTA
457/// where DELTA.GT.0 and D.GT.0
458
460{
461 Double_t t1 = 1./nsq;
462 Double_t t3=rsq-(radius*radius);
463 b = t1*rdotn;
464 Double_t c =t1*t3;
465 delta = b*b-c;
466 if (delta>0) {
467 delta=TMath::Sqrt(delta);
468 } else {
469 delta = -1;
470 }
471}
472
473////////////////////////////////////////////////////////////////////////////////
474/// Divide this tube shape belonging to volume "voldiv" into ndiv volumes
475/// called divname, from start position with the given step. Returns pointer
476/// to created division cell volume in case of Z divisions. For radial division
477/// creates all volumes with different shapes and returns pointer to volume that
478/// was divided. In case a wrong division axis is supplied, returns pointer to
479/// volume that was divided.
480
481TGeoVolume *TGeoTube::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
482 Double_t start, Double_t step)
483{
484 TGeoShape *shape; //--- shape to be created
485 TGeoVolume *vol; //--- division volume to be created
486 TGeoVolumeMulti *vmulti; //--- generic divided volume
487 TGeoPatternFinder *finder; //--- finder to be attached
488 TString opt = ""; //--- option to be attached
489 Int_t id;
490 Double_t end = start+ndiv*step;
491 switch (iaxis) {
492 case 1: //--- R division
493 finder = new TGeoPatternCylR(voldiv, ndiv, start, end);
494 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
495 voldiv->SetFinder(finder);
496 finder->SetDivIndex(voldiv->GetNdaughters());
497 for (id=0; id<ndiv; id++) {
498 shape = new TGeoTube(start+id*step, start+(id+1)*step, fDz);
499 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
500 vmulti->AddVolume(vol);
501 opt = "R";
502 voldiv->AddNodeOffset(vol, id, 0, opt.Data());
503 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
504 }
505 return vmulti;
506 case 2: //--- Phi division
507 finder = new TGeoPatternCylPhi(voldiv, ndiv, start, end);
508 voldiv->SetFinder(finder);
509 finder->SetDivIndex(voldiv->GetNdaughters());
510 shape = new TGeoTubeSeg(fRmin, fRmax, fDz, -step/2, step/2);
511 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
512 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
513 vmulti->AddVolume(vol);
514 opt = "Phi";
515 for (id=0; id<ndiv; id++) {
516 voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
517 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
518 }
519 return vmulti;
520 case 3: //--- Z division
521 finder = new TGeoPatternZ(voldiv, ndiv, start, start+ndiv*step);
522 voldiv->SetFinder(finder);
523 finder->SetDivIndex(voldiv->GetNdaughters());
524 shape = new TGeoTube(fRmin, fRmax, step/2);
525 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
526 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
527 vmulti->AddVolume(vol);
528 opt = "Z";
529 for (id=0; id<ndiv; id++) {
530 voldiv->AddNodeOffset(vol, id, start+step/2+id*step, opt.Data());
531 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
532 }
533 return vmulti;
534 default:
535 Error("Divide", "In shape %s wrong axis type for division", GetName());
536 return 0;
537 }
538}
539
540////////////////////////////////////////////////////////////////////////////////
541/// Returns name of axis IAXIS.
542
543const char *TGeoTube::GetAxisName(Int_t iaxis) const
544{
545 switch (iaxis) {
546 case 1:
547 return "R";
548 case 2:
549 return "PHI";
550 case 3:
551 return "Z";
552 default:
553 return "UNDEFINED";
554 }
555}
556
557////////////////////////////////////////////////////////////////////////////////
558/// Get range of shape for a given axis.
559
561{
562 xlo = 0;
563 xhi = 0;
564 Double_t dx = 0;
565 switch (iaxis) {
566 case 1:
567 xlo = fRmin;
568 xhi = fRmax;
569 dx = xhi-xlo;
570 return dx;
571 case 2:
572 xlo = 0;
573 xhi = 360;
574 dx = 360;
575 return dx;
576 case 3:
577 xlo = -fDz;
578 xhi = fDz;
579 dx = xhi-xlo;
580 return dx;
581 }
582 return dx;
583}
584
585////////////////////////////////////////////////////////////////////////////////
586/// Fill vector param[4] with the bounding cylinder parameters. The order
587/// is the following : Rmin, Rmax, Phi1, Phi2, dZ
588
590{
591 param[0] = fRmin; // Rmin
592 param[0] *= param[0];
593 param[1] = fRmax; // Rmax
594 param[1] *= param[1];
595 param[2] = 0.; // Phi1
596 param[3] = 360.; // Phi2
597}
598
599////////////////////////////////////////////////////////////////////////////////
600/// in case shape has some negative parameters, these has to be computed
601/// in order to fit the mother
602
604{
605 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
606 Double_t rmin, rmax, dz;
608 rmin = fRmin;
609 rmax = fRmax;
610 dz = fDz;
611 if (fDz<0) {
612 mother->GetAxisRange(3,xmin,xmax);
613 if (xmax<0) return 0;
614 dz=xmax;
615 }
616 mother->GetAxisRange(1,xmin,xmax);
617 if (fRmin<0) {
618 if (xmin<0) return 0;
619 rmin = xmin;
620 }
621 if (fRmax<0) {
622 if (xmax<=0) return 0;
623 rmax = xmax;
624 }
625
626 return (new TGeoTube(GetName(), rmin, rmax, dz));
627}
628
629////////////////////////////////////////////////////////////////////////////////
630/// print shape parameters
631
633{
634 printf("*** Shape %s: TGeoTube ***\n", GetName());
635 printf(" Rmin = %11.5f\n", fRmin);
636 printf(" Rmax = %11.5f\n", fRmax);
637 printf(" dz = %11.5f\n", fDz);
638 printf(" Bounding box:\n");
640}
641
642////////////////////////////////////////////////////////////////////////////////
643/// Creates a TBuffer3D describing *this* shape.
644/// Coordinates are in local reference frame.
645
647{
649 Int_t nbPnts = 4*n;
650 Int_t nbSegs = 8*n;
651 Int_t nbPols = 4*n;
652 if (!HasRmin()) {
653 nbPnts = 2*(n+1);
654 nbSegs = 5*n;
655 nbPols = 3*n;
656 }
658 nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols);
659 if (buff)
660 {
661 SetPoints(buff->fPnts);
662 SetSegsAndPols(*buff);
663 }
664
665 return buff;
666}
667
668////////////////////////////////////////////////////////////////////////////////
669/// Fill TBuffer3D structure for segments and polygons.
670
672{
673 Int_t i, j,indx;
675 Int_t c = (((buffer.fColor) %8) -1) * 4;
676 if (c < 0) c = 0;
677
678 if (HasRmin()) {
679 // circle segments:
680 // lower rmin circle: i=0, (0, n-1)
681 // lower rmax circle: i=1, (n, 2n-1)
682 // upper rmin circle: i=2, (2n, 3n-1)
683 // upper rmax circle: i=1, (3n, 4n-1)
684 for (i = 0; i < 4; i++) {
685 for (j = 0; j < n; j++) {
686 indx = 3*(i*n+j);
687 buffer.fSegs[indx ] = c;
688 buffer.fSegs[indx+1] = i*n+j;
689 buffer.fSegs[indx+2] = i*n+(j+1)%n;
690 }
691 }
692 // Z-parallel segments
693 // inner: i=4, (4n, 5n-1)
694 // outer: i=5, (5n, 6n-1)
695 for (i = 4; i < 6; i++) {
696 for (j = 0; j < n; j++) {
697 indx = 3*(i*n+j);
698 buffer.fSegs[indx ] = c+1;
699 buffer.fSegs[indx+1] = (i-4)*n+j;
700 buffer.fSegs[indx+2] = (i-2)*n+j;
701 }
702 }
703 // Radial segments
704 // lower: i=6, (6n, 7n-1)
705 // upper: i=7, (7n, 8n-1)
706 for (i = 6; i < 8; i++) {
707 for (j = 0; j < n; j++) {
708 indx = 3*(i*n+j);
709 buffer.fSegs[indx ] = c;
710 buffer.fSegs[indx+1] = 2*(i-6)*n+j;
711 buffer.fSegs[indx+2] = (2*(i-6)+1)*n+j;
712 }
713 }
714 // Polygons
715 i=0;
716 // Inner lateral (0, n-1)
717 for (j = 0; j < n; j++) {
718 indx = 6*(i*n+j);
719 buffer.fPols[indx ] = c;
720 buffer.fPols[indx+1] = 4;
721 buffer.fPols[indx+2] = j;
722 buffer.fPols[indx+3] = 4*n+(j+1)%n;
723 buffer.fPols[indx+4] = 2*n+j;
724 buffer.fPols[indx+5] = 4*n+j;
725 }
726 i=1;
727 // Outer lateral (n,2n-1)
728 for (j = 0; j < n; j++) {
729 indx = 6*(i*n+j);
730 buffer.fPols[indx ] = c+1;
731 buffer.fPols[indx+1] = 4;
732 buffer.fPols[indx+2] = n+j;
733 buffer.fPols[indx+3] = 5*n+j;
734 buffer.fPols[indx+4] = 3*n+j;
735 buffer.fPols[indx+5] = 5*n+(j+1)%n;
736 }
737 i=2;
738 // lower disc (2n, 3n-1)
739 for (j = 0; j < n; j++) {
740 indx = 6*(i*n+j);
741 buffer.fPols[indx ] = c;
742 buffer.fPols[indx+1] = 4;
743 buffer.fPols[indx+2] = j;
744 buffer.fPols[indx+3] = 6*n+j;
745 buffer.fPols[indx+4] = n+j;
746 buffer.fPols[indx+5] = 6*n+(j+1)%n;
747 }
748 i=3;
749 // upper disc (3n, 4n-1)
750 for (j = 0; j < n; j++) {
751 indx = 6*(i*n+j);
752 buffer.fPols[indx ] = c;
753 buffer.fPols[indx+1] = 4;
754 buffer.fPols[indx+2] = 2*n+j;
755 buffer.fPols[indx+3] = 7*n+(j+1)%n;
756 buffer.fPols[indx+4] = 3*n+j;
757 buffer.fPols[indx+5] = 7*n+j;
758 }
759 return;
760 }
761 // Rmin=0 tubes
762 // circle segments
763 // lower rmax circle: i=0, (0, n-1)
764 // upper rmax circle: i=1, (n, 2n-1)
765 for (i = 0; i < 2; i++) {
766 for (j = 0; j < n; j++) {
767 indx = 3*(i*n+j);
768 buffer.fSegs[indx ] = c;
769 buffer.fSegs[indx+1] = 2+i*n+j;
770 buffer.fSegs[indx+2] = 2+i*n+(j+1)%n;
771 }
772 }
773 // Z-parallel segments (2n,3n-1)
774 for (j = 0; j < n; j++) {
775 indx = 3*(2*n+j);
776 buffer.fSegs[indx ] = c+1;
777 buffer.fSegs[indx+1] = 2+j;
778 buffer.fSegs[indx+2] = 2+n+j;
779 }
780 // Radial segments
781 // Lower circle: i=3, (3n,4n-1)
782 // Upper circle: i=4, (4n,5n-1)
783 for (i = 3; i < 5; i++) {
784 for (j = 0; j < n; j++) {
785 indx = 3*(i*n+j);
786 buffer.fSegs[indx ] = c;
787 buffer.fSegs[indx+1] = i-3;
788 buffer.fSegs[indx+2] = 2+(i-3)*n+j;
789 }
790 }
791 // Polygons
792 // lateral (0,n-1)
793 for (j = 0; j < n; j++) {
794 indx = 6*j;
795 buffer.fPols[indx ] = c+1;
796 buffer.fPols[indx+1] = 4;
797 buffer.fPols[indx+2] = j;
798 buffer.fPols[indx+3] = 2*n+j;
799 buffer.fPols[indx+4] = n+j;
800 buffer.fPols[indx+5] = 2*n+(j+1)%n;
801 }
802 // bottom triangles (n,2n-1)
803 for (j = 0; j < n; j++) {
804 indx = 6*n + 5*j;
805 buffer.fPols[indx ] = c;
806 buffer.fPols[indx+1] = 3;
807 buffer.fPols[indx+2] = j;
808 buffer.fPols[indx+3] = 3*n+(j+1)%n;
809 buffer.fPols[indx+4] = 3*n+j;
810 }
811 // top triangles (2n,3n-1)
812 for (j = 0; j < n; j++) {
813 indx = 6*n + 5*n + 5*j;
814 buffer.fPols[indx ] = c;
815 buffer.fPols[indx+1] = 3;
816 buffer.fPols[indx+2] = n+j;
817 buffer.fPols[indx+3] = 4*n+j;
818 buffer.fPols[indx+4] = 4*n+(j+1)%n;
819 }
820}
821
822////////////////////////////////////////////////////////////////////////////////
823/// computes the closest distance from given point to this shape, according
824/// to option. The matching point on the shape is stored in spoint.
825
827{
828#ifndef NEVER
829 Double_t r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
830 Double_t safe, safrmin, safrmax;
831 if (in) {
832 safe = fDz-TMath::Abs(point[2]); // positive if inside
833 if (fRmin>1E-10) {
834 safrmin = r-fRmin;
835 if (safrmin < safe) safe = safrmin;
836 }
837 safrmax = fRmax-r;
838 if (safrmax < safe) safe = safrmax;
839 } else {
840 safe = -fDz+TMath::Abs(point[2]);
841 if (fRmin>1E-10) {
842 safrmin = -r+fRmin;
843 if (safrmin > safe) safe = safrmin;
844 }
845 safrmax = -fRmax+r;
846 if (safrmax > safe) safe = safrmax;
847 }
848 return safe;
849#else
850 Double_t saf[3];
851 Double_t rsq = point[0]*point[0]+point[1]*point[1];
852 Double_t r = TMath::Sqrt(rsq);
853 saf[0] = fDz-TMath::Abs(point[2]); // positive if inside
854 saf[1] = (fRmin>1E-10)?(r-fRmin):TGeoShape::Big();
855 saf[2] = fRmax-r;
856 if (in) return saf[TMath::LocMin(3,saf)];
857 for (Int_t i=0; i<3; i++) saf[i]=-saf[i];
858 return saf[TMath::LocMax(3,saf)];
859#endif
860}
861
862////////////////////////////////////////////////////////////////////////////////
863/// computes the closest distance from given point to this shape, according
864/// to option. The matching point on the shape is stored in spoint.
865
866Double_t TGeoTube::SafetyS(const Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Int_t skipz)
867{
868 Double_t saf[3];
869 Double_t rsq = point[0]*point[0]+point[1]*point[1];
870 Double_t r = TMath::Sqrt(rsq);
871 switch (skipz) {
872 case 1: // skip lower Z plane
873 saf[0] = dz - point[2];
874 break;
875 case 2: // skip upper Z plane
876 saf[0] = dz + point[2];
877 break;
878 case 3: // skip both
879 saf[0] = TGeoShape::Big();
880 break;
881 default:
882 saf[0] = dz-TMath::Abs(point[2]);
883 }
884 saf[1] = (rmin>1E-10)?(r-rmin):TGeoShape::Big();
885 saf[2] = rmax-r;
886// printf("saf0=%g saf1=%g saf2=%g in=%d skipz=%d\n", saf[0],saf[1],saf[2],in,skipz);
887 if (in) return saf[TMath::LocMin(3,saf)];
888 for (Int_t i=0; i<3; i++) saf[i]=-saf[i];
889 return saf[TMath::LocMax(3,saf)];
890}
891
892////////////////////////////////////////////////////////////////////////////////
893/// Save a primitive as a C++ statement(s) on output stream "out".
894
895void TGeoTube::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
896{
898 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
899 out << " rmin = " << fRmin << ";" << std::endl;
900 out << " rmax = " << fRmax << ";" << std::endl;
901 out << " dz = " << fDz << ";" << std::endl;
902 out << " TGeoShape *" << GetPointerName() << " = new TGeoTube(\"" << GetName() << "\",rmin,rmax,dz);" << std::endl;
904}
905
906////////////////////////////////////////////////////////////////////////////////
907/// Set tube dimensions.
908
910{
911 fRmin = rmin;
912 fRmax = rmax;
913 fDz = dz;
914 if (fRmin>0 && fRmax>0 && fRmin>=fRmax)
915 Error("SetTubeDimensions", "In shape %s wrong rmin=%g rmax=%g", GetName(), rmin,rmax);
916}
917
918////////////////////////////////////////////////////////////////////////////////
919/// Set tube dimensions starting from a list.
920
922{
923 Double_t rmin = param[0];
924 Double_t rmax = param[1];
925 Double_t dz = param[2];
926 SetTubeDimensions(rmin, rmax, dz);
927}
928
929////////////////////////////////////////////////////////////////////////////////
930/// Fills array with n random points located on the line segments of the shape mesh.
931/// The output array must be provided with a length of minimum 3*npoints. Returns
932/// true if operation is implemented.
933
935{
936 if (npoints > (npoints/2)*2) {
937 Error("GetPointsOnSegments","Npoints must be even number");
938 return kFALSE;
939 }
940 Int_t nc = 0;
941 if (HasRmin()) nc = (Int_t)TMath::Sqrt(0.5*npoints);
942 else nc = (Int_t)TMath::Sqrt(1.*npoints);
943 Double_t dphi = TMath::TwoPi()/nc;
944 Double_t phi = 0;
945 Int_t ntop = 0;
946 if (HasRmin()) ntop = npoints/2 - nc*(nc-1);
947 else ntop = npoints - nc*(nc-1);
948 Double_t dz = 2*fDz/(nc-1);
949 Double_t z = 0;
950 Int_t icrt = 0;
951 Int_t nphi = nc;
952 // loop z sections
953 for (Int_t i=0; i<nc; i++) {
954 if (i == (nc-1)) nphi = ntop;
955 z = -fDz + i*dz;
956 // loop points on circle sections
957 for (Int_t j=0; j<nphi; j++) {
958 phi = j*dphi;
959 if (HasRmin()) {
960 array[icrt++] = fRmin * TMath::Cos(phi);
961 array[icrt++] = fRmin * TMath::Sin(phi);
962 array[icrt++] = z;
963 }
964 array[icrt++] = fRmax * TMath::Cos(phi);
965 array[icrt++] = fRmax * TMath::Sin(phi);
966 array[icrt++] = z;
967 }
968 }
969 return kTRUE;
970}
971
972////////////////////////////////////////////////////////////////////////////////
973/// create tube mesh points
974
976{
977 Double_t dz;
978 Int_t j, n;
980 Double_t dphi = 360./n;
981 Double_t phi = 0;
982 dz = fDz;
983 Int_t indx = 0;
984 if (points) {
985 if (HasRmin()) {
986 // 4*n points
987 // (0,n-1) lower rmin circle
988 // (2n, 3n-1) upper rmin circle
989 for (j = 0; j < n; j++) {
990 phi = j*dphi*TMath::DegToRad();
991 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
992 indx++;
993 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
994 indx++;
995 points[indx+6*n] = dz;
996 points[indx] =-dz;
997 indx++;
998 }
999 // (n, 2n-1) lower rmax circle
1000 // (3n, 4n-1) upper rmax circle
1001 for (j = 0; j < n; j++) {
1002 phi = j*dphi*TMath::DegToRad();
1003 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
1004 indx++;
1005 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
1006 indx++;
1007 points[indx+6*n]= dz;
1008 points[indx] =-dz;
1009 indx++;
1010 }
1011 } else {
1012 // centers of lower/upper circles (0,1)
1013 points[indx++] = 0.;
1014 points[indx++] = 0.;
1015 points[indx++] = -dz;
1016 points[indx++] = 0.;
1017 points[indx++] = 0.;
1018 points[indx++] = dz;
1019 // lower rmax circle (2, 2+n-1)
1020 // upper rmax circle (2+n, 2+2n-1)
1021 for (j = 0; j < n; j++) {
1022 phi = j*dphi*TMath::DegToRad();
1023 points[indx+3*n] = points[indx] = fRmax * TMath::Cos(phi);
1024 indx++;
1025 points[indx+3*n] = points[indx] = fRmax * TMath::Sin(phi);
1026 indx++;
1027 points[indx+3*n]= dz;
1028 points[indx] =-dz;
1029 indx++;
1030 }
1031 }
1032 }
1033}
1034
1035////////////////////////////////////////////////////////////////////////////////
1036/// create tube mesh points
1037
1039{
1040 Double_t dz;
1041 Int_t j, n;
1043 Double_t dphi = 360./n;
1044 Double_t phi = 0;
1045 dz = fDz;
1046 Int_t indx = 0;
1047 if (points) {
1048 if (HasRmin()) {
1049 // 4*n points
1050 // (0,n-1) lower rmin circle
1051 // (2n, 3n-1) upper rmin circle
1052 for (j = 0; j < n; j++) {
1053 phi = j*dphi*TMath::DegToRad();
1054 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
1055 indx++;
1056 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
1057 indx++;
1058 points[indx+6*n] = dz;
1059 points[indx] =-dz;
1060 indx++;
1061 }
1062 // (n, 2n-1) lower rmax circle
1063 // (3n, 4n-1) upper rmax circle
1064 for (j = 0; j < n; j++) {
1065 phi = j*dphi*TMath::DegToRad();
1066 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
1067 indx++;
1068 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
1069 indx++;
1070 points[indx+6*n]= dz;
1071 points[indx] =-dz;
1072 indx++;
1073 }
1074 } else {
1075 // centers of lower/upper circles (0,1)
1076 points[indx++] = 0.;
1077 points[indx++] = 0.;
1078 points[indx++] = -dz;
1079 points[indx++] = 0.;
1080 points[indx++] = 0.;
1081 points[indx++] = dz;
1082 // lower rmax circle (2, 2+n-1)
1083 // upper rmax circle (2+n, 2+2n-1)
1084 for (j = 0; j < n; j++) {
1085 phi = j*dphi*TMath::DegToRad();
1086 points[indx+3*n] = points[indx] = fRmax * TMath::Cos(phi);
1087 indx++;
1088 points[indx+3*n] = points[indx] = fRmax * TMath::Sin(phi);
1089 indx++;
1090 points[indx+3*n]= dz;
1091 points[indx] =-dz;
1092 indx++;
1093 }
1094 }
1095 }
1096}
1097
1098////////////////////////////////////////////////////////////////////////////////
1099/// Return number of vertices of the mesh representation
1100
1102{
1104 Int_t numPoints = n*4;
1105 if (!HasRmin()) numPoints = 2*(n+1);
1106 return numPoints;
1107}
1108
1109////////////////////////////////////////////////////////////////////////////////
1110/// Returns numbers of vertices, segments and polygons composing the shape mesh.
1111
1112void TGeoTube::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
1113{
1115 nvert = n*4;
1116 nsegs = n*8;
1117 npols = n*4;
1118 if (!HasRmin()) {
1119 nvert = 2*(n+1);
1120 nsegs = 5*n;
1121 npols = 3*n;
1122 } else {
1123 nvert = n*4;
1124 nsegs = n*8;
1125 npols = n*4;
1126 }
1127}
1128
1129////////////////////////////////////////////////////////////////////////////////
1130/// fill size of this 3-D object
1131
1133{
1134}
1135
1136////////////////////////////////////////////////////////////////////////////////
1137/// Fills a static 3D buffer and returns a reference.
1138
1139const TBuffer3D & TGeoTube::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
1140{
1141 static TBuffer3DTube buffer;
1142 TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
1143
1144 if (reqSections & TBuffer3D::kShapeSpecific) {
1145 buffer.fRadiusInner = fRmin;
1146 buffer.fRadiusOuter = fRmax;
1147 buffer.fHalfLength = fDz;
1149 }
1150 if (reqSections & TBuffer3D::kRawSizes) {
1152 Int_t nbPnts = 4*n;
1153 Int_t nbSegs = 8*n;
1154 Int_t nbPols = 4*n;
1155 if (!HasRmin()) {
1156 nbPnts = 2*(n+1);
1157 nbSegs = 5*n;
1158 nbPols = 3*n;
1159 }
1160 if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
1162 }
1163 }
1164 if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
1165 SetPoints(buffer.fPnts);
1166 if (!buffer.fLocalFrame) {
1167 TransformPoints(buffer.fPnts, buffer.NbPnts());
1168 }
1169 SetSegsAndPols(buffer);
1171 }
1172
1173 return buffer;
1174}
1175
1176////////////////////////////////////////////////////////////////////////////////
1177/// Check the inside status for each of the points in the array.
1178/// Input: Array of point coordinates + vector size
1179/// Output: Array of Booleans for the inside of each point
1180
1181void TGeoTube::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
1182{
1183 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
1184}
1185
1186////////////////////////////////////////////////////////////////////////////////
1187/// Compute the normal for an array o points so that norm.dot.dir is positive
1188/// Input: Arrays of point coordinates and directions + vector size
1189/// Output: Array of normal directions
1190
1191void TGeoTube::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
1192{
1193 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
1194}
1195
1196////////////////////////////////////////////////////////////////////////////////
1197/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1198
1199void TGeoTube::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1200{
1201 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
1202}
1203
1204////////////////////////////////////////////////////////////////////////////////
1205/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1206
1207void TGeoTube::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1208{
1209 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
1210}
1211
1212////////////////////////////////////////////////////////////////////////////////
1213/// Compute safe distance from each of the points in the input array.
1214/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1215/// Output: Safety values
1216
1217void TGeoTube::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
1218{
1219 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
1220}
1221
1223
1224////////////////////////////////////////////////////////////////////////////////
1225/// Default constructor
1226
1228 :TGeoTube(),
1229 fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1230{
1232}
1233
1234////////////////////////////////////////////////////////////////////////////////
1235/// Default constructor specifying minimum and maximum radius.
1236/// The segment will be from phiStart to phiEnd expressed in degree.
1237
1239 Double_t phiStart, Double_t phiEnd)
1240 :TGeoTube(rmin, rmax, dz),
1241 fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1242{
1244 SetTubsDimensions(rmin, rmax, dz, phiStart, phiEnd);
1245 ComputeBBox();
1246}
1247
1248////////////////////////////////////////////////////////////////////////////////
1249/// Default constructor specifying minimum and maximum radius
1250/// The segment will be from phiStart to phiEnd expressed in degree.
1251
1253 Double_t phiStart, Double_t phiEnd)
1254 :TGeoTube(name, rmin, rmax, dz)
1255{
1257 SetTubsDimensions(rmin, rmax, dz, phiStart, phiEnd);
1258 ComputeBBox();
1259}
1260
1261////////////////////////////////////////////////////////////////////////////////
1262/// Default constructor specifying minimum and maximum radius
1263/// - param[0] = Rmin
1264/// - param[1] = Rmax
1265/// - param[2] = dz
1266/// - param[3] = phi1
1267/// - param[4] = phi2
1268
1270 :TGeoTube(0, 0, 0)
1271{
1273 SetDimensions(param);
1274 ComputeBBox();
1275}
1276
1277////////////////////////////////////////////////////////////////////////////////
1278/// destructor
1279
1281{
1282}
1283
1284////////////////////////////////////////////////////////////////////////////////
1285/// Function called after streaming an object of this class.
1286
1288{
1290}
1291
1292////////////////////////////////////////////////////////////////////////////////
1293/// Init frequently used trigonometric values
1294
1296{
1299 fC1 = TMath::Cos(phi1);
1300 fS1 = TMath::Sin(phi1);
1301 fC2 = TMath::Cos(phi2);
1302 fS2 = TMath::Sin(phi2);
1303 Double_t fio = 0.5*(phi1+phi2);
1304 fCm = TMath::Cos(fio);
1305 fSm = TMath::Sin(fio);
1306 Double_t dfi = 0.5*(phi2-phi1);
1307 fCdfi = TMath::Cos(dfi);
1308}
1309
1310////////////////////////////////////////////////////////////////////////////////
1311/// Computes capacity of the shape in [length^3]
1312
1314{
1316}
1317
1318////////////////////////////////////////////////////////////////////////////////
1319/// Computes capacity of the shape in [length^3]
1320
1322{
1323 Double_t capacity = TMath::Abs(phiEnd-phiStart)*TMath::DegToRad()*(rmax*rmax-rmin*rmin)*dz;
1324 return capacity;
1325}
1326
1327////////////////////////////////////////////////////////////////////////////////
1328/// compute bounding box of the tube segment
1329
1331{
1332 Double_t xc[4];
1333 Double_t yc[4];
1334 xc[0] = fRmax*fC1;
1335 yc[0] = fRmax*fS1;
1336 xc[1] = fRmax*fC2;
1337 yc[1] = fRmax*fS2;
1338 xc[2] = fRmin*fC1;
1339 yc[2] = fRmin*fS1;
1340 xc[3] = fRmin*fC2;
1341 yc[3] = fRmin*fS2;
1342
1343 Double_t xmin = xc[TMath::LocMin(4, &xc[0])];
1344 Double_t xmax = xc[TMath::LocMax(4, &xc[0])];
1345 Double_t ymin = yc[TMath::LocMin(4, &yc[0])];
1346 Double_t ymax = yc[TMath::LocMax(4, &yc[0])];
1347
1348 Double_t dp = fPhi2-fPhi1;
1349 if (dp<0) dp+=360;
1350 Double_t ddp = -fPhi1;
1351 if (ddp<0) ddp+= 360;
1352 if (ddp>360) ddp-=360;
1353 if (ddp<=dp) xmax = fRmax;
1354 ddp = 90-fPhi1;
1355 if (ddp<0) ddp+= 360;
1356 if (ddp>360) ddp-=360;
1357 if (ddp<=dp) ymax = fRmax;
1358 ddp = 180-fPhi1;
1359 if (ddp<0) ddp+= 360;
1360 if (ddp>360) ddp-=360;
1361 if (ddp<=dp) xmin = -fRmax;
1362 ddp = 270-fPhi1;
1363 if (ddp<0) ddp+= 360;
1364 if (ddp>360) ddp-=360;
1365 if (ddp<=dp) ymin = -fRmax;
1366 fOrigin[0] = (xmax+xmin)/2;
1367 fOrigin[1] = (ymax+ymin)/2;
1368 fOrigin[2] = 0;
1369 fDX = (xmax-xmin)/2;
1370 fDY = (ymax-ymin)/2;
1371 fDZ = fDz;
1372}
1373
1374////////////////////////////////////////////////////////////////////////////////
1375/// Compute normal to closest surface from POINT.
1376
1377void TGeoTubeSeg::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
1378{
1379 Double_t saf[3];
1380 Double_t rsq = point[0]*point[0]+point[1]*point[1];
1381 Double_t r = TMath::Sqrt(rsq);
1382 saf[0] = TMath::Abs(fDz-TMath::Abs(point[2]));
1383 saf[1] = (fRmin>1E-10)?TMath::Abs(r-fRmin):TGeoShape::Big();
1384 saf[2] = TMath::Abs(fRmax-r);
1385 Int_t i = TMath::LocMin(3,saf);
1386 if (((fPhi2-fPhi1)<360.) && TGeoShape::IsCloseToPhi(saf[i], point,fC1,fS1,fC2,fS2)) {
1387 TGeoShape::NormalPhi(point,dir,norm,fC1,fS1,fC2,fS2);
1388 return;
1389 }
1390 if (i==0) {
1391 norm[0] = norm[1] = 0.;
1392 norm[2] = TMath::Sign(1.,dir[2]);
1393 return;
1394 };
1395 norm[2] = 0;
1396 Double_t phi = TMath::ATan2(point[1], point[0]);
1397 norm[0] = TMath::Cos(phi);
1398 norm[1] = TMath::Sin(phi);
1399 if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
1400 norm[0] = -norm[0];
1401 norm[1] = -norm[1];
1402 }
1403}
1404
1405////////////////////////////////////////////////////////////////////////////////
1406/// Compute normal to closest surface from POINT.
1407
1408void TGeoTubeSeg::ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm,
1409 Double_t rmin, Double_t rmax, Double_t /*dz*/,
1411{
1412 Double_t saf[2];
1413 Double_t rsq = point[0]*point[0]+point[1]*point[1];
1414 Double_t r = TMath::Sqrt(rsq);
1415 saf[0] = (rmin>1E-10)?TMath::Abs(r-rmin):TGeoShape::Big();
1416 saf[1] = TMath::Abs(rmax-r);
1417 Int_t i = TMath::LocMin(2,saf);
1418 if (TGeoShape::IsCloseToPhi(saf[i], point,c1,s1,c2,s2)) {
1419 TGeoShape::NormalPhi(point,dir,norm,c1,s1,c2,s2);
1420 return;
1421 }
1422 norm[2] = 0;
1423 Double_t phi = TMath::ATan2(point[1], point[0]);
1424 norm[0] = TMath::Cos(phi);
1425 norm[1] = TMath::Sin(phi);
1426 if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
1427 norm[0] = -norm[0];
1428 norm[1] = -norm[1];
1429 }
1430}
1431
1432////////////////////////////////////////////////////////////////////////////////
1433/// test if point is inside this tube segment
1434/// first check if point is inside the tube
1435
1437{
1438 if (!TGeoTube::Contains(point)) return kFALSE;
1439 return IsInPhiRange(point, fPhi1, fPhi2);
1440}
1441
1442////////////////////////////////////////////////////////////////////////////////
1443/// compute closest distance from point px,py to each corner
1444
1446{
1448 const Int_t numPoints = 4*n;
1449 return ShapeDistancetoPrimitive(numPoints, px, py);
1450}
1451
1452////////////////////////////////////////////////////////////////////////////////
1453/// Compute distance from inside point to surface of the tube segment (static)
1454/// Boundary safe algorithm.
1455/// Do Z
1456
1459{
1460 Double_t stube = TGeoTube::DistFromInsideS(point,dir,rmin,rmax,dz);
1461 if (stube<=0) return 0.0;
1462 Double_t sfmin = TGeoShape::Big();
1463 Double_t rsq = point[0]*point[0]+point[1]*point[1];
1464 Double_t r = TMath::Sqrt(rsq);
1465 Double_t cpsi=point[0]*cm+point[1]*sm;
1466 if (cpsi>r*cdfi+TGeoShape::Tolerance()) {
1467 sfmin = TGeoShape::DistToPhiMin(point, dir, s1, c1, s2, c2, sm, cm);
1468 return TMath::Min(stube,sfmin);
1469 }
1470 // Point on the phi boundary or outside
1471 // which one: phi1 or phi2
1472 Double_t ddotn, xi, yi;
1473 if (TMath::Abs(point[1]-s1*r) < TMath::Abs(point[1]-s2*r)) {
1474 ddotn = s1*dir[0]-c1*dir[1];
1475 if (ddotn>=0) return 0.0;
1476 ddotn = -s2*dir[0]+c2*dir[1];
1477 if (ddotn<=0) return stube;
1478 sfmin = s2*point[0]-c2*point[1];
1479 if (sfmin<=0) return stube;
1480 sfmin /= ddotn;
1481 if (sfmin >= stube) return stube;
1482 xi = point[0]+sfmin*dir[0];
1483 yi = point[1]+sfmin*dir[1];
1484 if (yi*cm-xi*sm<0) return stube;
1485 return sfmin;
1486 }
1487 ddotn = -s2*dir[0]+c2*dir[1];
1488 if (ddotn>=0) return 0.0;
1489 ddotn = s1*dir[0]-c1*dir[1];
1490 if (ddotn<=0) return stube;
1491 sfmin = -s1*point[0]+c1*point[1];
1492 if (sfmin<=0) return stube;
1493 sfmin /= ddotn;
1494 if (sfmin >= stube) return stube;
1495 xi = point[0]+sfmin*dir[0];
1496 yi = point[1]+sfmin*dir[1];
1497 if (yi*cm-xi*sm>0) return stube;
1498 return sfmin;
1499}
1500
1501////////////////////////////////////////////////////////////////////////////////
1502/// Compute distance from inside point to surface of the tube segment
1503/// Boundary safe algorithm.
1504
1505Double_t TGeoTubeSeg::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1506{
1507 if (iact<3 && safe) {
1508 *safe = SafetyS(point, kTRUE, fRmin, fRmax, fDz, fPhi1, fPhi2);
1509 if (iact==0) return TGeoShape::Big();
1510 if ((iact==1) && (*safe>step)) return TGeoShape::Big();
1511 }
1512 if ((fPhi2-fPhi1)>=360.) return TGeoTube::DistFromInsideS(point,dir,fRmin,fRmax,fDz);
1513
1514 // compute distance to surface
1516}
1517
1518////////////////////////////////////////////////////////////////////////////////
1519/// Static method to compute distance to arbitrary tube segment from outside point
1520/// Boundary safe algorithm.
1521
1524 Double_t cm, Double_t sm, Double_t cdfi)
1525{
1526 Double_t r2, cpsi;
1527 // check Z planes
1528 Double_t xi, yi, zi;
1529 zi = dz - TMath::Abs(point[2]);
1530 Double_t rmaxsq = rmax*rmax;
1531 Double_t rminsq = rmin*rmin;
1533 Double_t snxt=TGeoShape::Big();
1534 Bool_t in = kFALSE;
1535 Bool_t inz = (zi<0)?kFALSE:kTRUE;
1536 if (!inz) {
1537 if (point[2]*dir[2]>=0) return TGeoShape::Big();
1538 s = -zi/TMath::Abs(dir[2]);
1539 xi = point[0]+s*dir[0];
1540 yi = point[1]+s*dir[1];
1541 r2=xi*xi+yi*yi;
1542 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1543 cpsi=(xi*cm+yi*sm)/TMath::Sqrt(r2);
1544 if (cpsi>=cdfi) return s;
1545 }
1546 }
1547
1548 // check outer cyl. surface
1549 Double_t rsq = point[0]*point[0]+point[1]*point[1];
1550 Double_t r = TMath::Sqrt(rsq);
1551 Double_t nsq=dir[0]*dir[0]+dir[1]*dir[1];
1552 Double_t rdotn=point[0]*dir[0]+point[1]*dir[1];
1553 Double_t b,d;
1554 Bool_t inrmax = kFALSE;
1555 Bool_t inrmin = kFALSE;
1556 Bool_t inphi = kFALSE;
1557 if (rsq<=rmaxsq+TGeoShape::Tolerance()) inrmax = kTRUE;
1558 if (rsq>=rminsq-TGeoShape::Tolerance()) inrmin = kTRUE;
1559 cpsi=point[0]*cm+point[1]*sm;
1560 if (cpsi>r*cdfi-TGeoShape::Tolerance()) inphi = kTRUE;
1561 in = inz & inrmin & inrmax & inphi;
1562 // If inside, we are most likely on a boundary within machine precision.
1563 if (in) {
1564 Bool_t checkout = kFALSE;
1565 Double_t safphi = (cpsi-r*cdfi)*TMath::Sqrt(1.-cdfi*cdfi);
1566// Double_t sch, cch;
1567 // check if on Z boundaries
1568 if (zi<rmax-r) {
1569 if (TGeoShape::IsSameWithinTolerance(rmin,0) || (zi<r-rmin)) {
1570 if (zi<safphi) {
1571 if (point[2]*dir[2]<0) return 0.0;
1572 return TGeoShape::Big();
1573 }
1574 }
1575 }
1576 if ((rmaxsq-rsq) < (rsq-rminsq)) checkout = kTRUE;
1577 // check if on Rmax boundary
1578 if (checkout && (rmax-r<safphi)) {
1579 if (rdotn>=0) return TGeoShape::Big();
1580 return 0.0;
1581 }
1582 if (TMath::Abs(nsq)<TGeoShape::Tolerance()) return TGeoShape::Big();
1583 // check if on phi boundary
1584 if (TGeoShape::IsSameWithinTolerance(rmin,0) || (safphi<r-rmin)) {
1585 // We may cross again a phi of rmin boundary
1586 // check first if we are on phi1 or phi2
1587 Double_t un;
1588 if (point[0]*c1 + point[1]*s1 > point[0]*c2 + point[1]*s2) {
1589 un = dir[0]*s1-dir[1]*c1;
1590 if (un < 0) return 0.0;
1591 if (cdfi>=0) return TGeoShape::Big();
1592 un = -dir[0]*s2+dir[1]*c2;
1593 if (un<0) {
1594 s = -point[0]*s2+point[1]*c2;
1595 if (s>0) {
1596 s /= (-un);
1597 zi = point[2]+s*dir[2];
1598 if (TMath::Abs(zi)<=dz) {
1599 xi = point[0]+s*dir[0];
1600 yi = point[1]+s*dir[1];
1601 r2=xi*xi+yi*yi;
1602 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1603 if ((yi*cm-xi*sm)>0) return s;
1604 }
1605 }
1606 }
1607 }
1608 } else {
1609 un = -dir[0]*s2+dir[1]*c2;
1610 if (un < 0) return 0.0;
1611 if (cdfi>=0) return TGeoShape::Big();
1612 un = dir[0]*s1-dir[1]*c1;
1613 if (un<0) {
1614 s = point[0]*s1-point[1]*c1;
1615 if (s>0) {
1616 s /= (-un);
1617 zi = point[2]+s*dir[2];
1618 if (TMath::Abs(zi)<=dz) {
1619 xi = point[0]+s*dir[0];
1620 yi = point[1]+s*dir[1];
1621 r2=xi*xi+yi*yi;
1622 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1623 if ((yi*cm-xi*sm)<0) return s;
1624 }
1625 }
1626 }
1627 }
1628 }
1629 // We may also cross rmin, (+) solution
1630 if (rdotn>=0) return TGeoShape::Big();
1631 if (cdfi>=0) return TGeoShape::Big();
1632 DistToTube(rsq, nsq, rdotn, rmin, b, d);
1633 if (d>0) {
1634 s=-b+d;
1635 if (s>0) {
1636 zi=point[2]+s*dir[2];
1637 if (TMath::Abs(zi)<=dz) {
1638 xi=point[0]+s*dir[0];
1639 yi=point[1]+s*dir[1];
1640 if ((xi*cm+yi*sm) >= rmin*cdfi) return s;
1641 }
1642 }
1643 }
1644 return TGeoShape::Big();
1645 }
1646 // we are on rmin boundary: we may cross again rmin or a phi facette
1647 if (rdotn>=0) return 0.0;
1648 DistToTube(rsq, nsq, rdotn, rmin, b, d);
1649 if (d>0) {
1650 s=-b+d;
1651 if (s>0) {
1652 zi=point[2]+s*dir[2];
1653 if (TMath::Abs(zi)<=dz) {
1654 // now check phi range
1655 xi=point[0]+s*dir[0];
1656 yi=point[1]+s*dir[1];
1657 r2=xi*xi+yi*yi;
1658 if ((xi*cm+yi*sm) >= rmin*cdfi) return s;
1659 // now we really have to check any phi crossing
1660 Double_t un=-dir[0]*s1+dir[1]*c1;
1661 if (un > 0) {
1662 s=point[0]*s1-point[1]*c1;
1663 if (s>=0) {
1664 s /= un;
1665 zi=point[2]+s*dir[2];
1666 if (TMath::Abs(zi)<=dz) {
1667 xi=point[0]+s*dir[0];
1668 yi=point[1]+s*dir[1];
1669 r2=xi*xi+yi*yi;
1670 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1671 if ((yi*cm-xi*sm)<=0) {
1672 if (s<snxt) snxt=s;
1673 }
1674 }
1675 }
1676 }
1677 }
1678 un=dir[0]*s2-dir[1]*c2;
1679 if (un > 0) {
1680 s=(point[1]*c2-point[0]*s2)/un;
1681 if (s>=0 && s<snxt) {
1682 zi=point[2]+s*dir[2];
1683 if (TMath::Abs(zi)<=dz) {
1684 xi=point[0]+s*dir[0];
1685 yi=point[1]+s*dir[1];
1686 r2=xi*xi+yi*yi;
1687 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1688 if ((yi*cm-xi*sm)>=0) {
1689 return s;
1690 }
1691 }
1692 }
1693 }
1694 }
1695 return snxt;
1696 }
1697 }
1698 }
1699 return TGeoShape::Big();
1700 }
1701 // only r>rmax has to be considered
1702 if (TMath::Abs(nsq)<TGeoShape::Tolerance()) return TGeoShape::Big();
1703 if (rsq>=rmax*rmax) {
1704 if (rdotn>=0) return TGeoShape::Big();
1705 TGeoTube::DistToTube(rsq, nsq, rdotn, rmax, b, d);
1706 if (d>0) {
1707 s=-b-d;
1708 if (s>0) {
1709 zi=point[2]+s*dir[2];
1710 if (TMath::Abs(zi)<=dz) {
1711 xi=point[0]+s*dir[0];
1712 yi=point[1]+s*dir[1];
1713 cpsi = xi*cm+yi*sm;
1714 if (cpsi>=rmax*cdfi) return s;
1715 }
1716 }
1717 }
1718 }
1719 // check inner cylinder
1720 if (rmin>0) {
1721 TGeoTube::DistToTube(rsq, nsq, rdotn, rmin, b, d);
1722 if (d>0) {
1723 s=-b+d;
1724 if (s>0) {
1725 zi=point[2]+s*dir[2];
1726 if (TMath::Abs(zi)<=dz) {
1727 xi=point[0]+s*dir[0];
1728 yi=point[1]+s*dir[1];
1729 cpsi = xi*cm+yi*sm;
1730 if (cpsi>=rmin*cdfi) snxt=s;
1731 }
1732 }
1733 }
1734 }
1735 // check phi planes
1736 Double_t un=-dir[0]*s1+dir[1]*c1;
1737 if (un > 0) {
1738 s=point[0]*s1-point[1]*c1;
1739 if (s>=0) {
1740 s /= un;
1741 zi=point[2]+s*dir[2];
1742 if (TMath::Abs(zi)<=dz) {
1743 xi=point[0]+s*dir[0];
1744 yi=point[1]+s*dir[1];
1745 r2=xi*xi+yi*yi;
1746 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1747 if ((yi*cm-xi*sm)<=0) {
1748 if (s<snxt) snxt=s;
1749 }
1750 }
1751 }
1752 }
1753 }
1754 un=dir[0]*s2-dir[1]*c2;
1755 if (un > 0) {
1756 s=point[1]*c2-point[0]*s2;
1757 if (s>=0) {
1758 s /= un;
1759 zi=point[2]+s*dir[2];
1760 if (TMath::Abs(zi)<=dz) {
1761 xi=point[0]+s*dir[0];
1762 yi=point[1]+s*dir[1];
1763 r2=xi*xi+yi*yi;
1764 if ((rminsq<=r2) && (r2<=rmaxsq)) {
1765 if ((yi*cm-xi*sm)>=0) {
1766 if (s<snxt) snxt=s;
1767 }
1768 }
1769 }
1770 }
1771 }
1772 return snxt;
1773}
1774
1775////////////////////////////////////////////////////////////////////////////////
1776/// compute distance from outside point to surface of the tube segment
1777/// fist localize point w.r.t tube
1778
1779Double_t TGeoTubeSeg::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1780{
1781 if (iact<3 && safe) {
1782 *safe = SafetyS(point, kFALSE, fRmin, fRmax, fDz, fPhi1, fPhi2);
1783 if (iact==0) return TGeoShape::Big();
1784 if ((iact==1) && (step<=*safe)) return TGeoShape::Big();
1785 }
1786// Check if the bounding box is crossed within the requested distance
1787 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
1788 if (sdist>=step) return TGeoShape::Big();
1789 if ((fPhi2-fPhi1)>=360.) return TGeoTube::DistFromOutsideS(point,dir,fRmin,fRmax,fDz);
1790
1791 // find distance to shape
1792 return TGeoTubeSeg::DistFromOutsideS(point, dir, fRmin, fRmax, fDz, fC1, fS1, fC2, fS2, fCm, fSm, fCdfi);
1793}
1794
1795////////////////////////////////////////////////////////////////////////////////
1796/// Divide this tube segment shape belonging to volume "voldiv" into ndiv volumes
1797/// called divname, from start position with the given step. Returns pointer
1798/// to created division cell volume in case of Z divisions. For radialdivision
1799/// creates all volumes with different shapes and returns pointer to volume that
1800/// was divided. In case a wrong division axis is supplied, returns pointer to
1801/// volume that was divided.
1802
1803TGeoVolume *TGeoTubeSeg::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
1804 Double_t start, Double_t step)
1805{
1806 TGeoShape *shape; //--- shape to be created
1807 TGeoVolume *vol; //--- division volume to be created
1808 TGeoVolumeMulti *vmulti; //--- generic divided volume
1809 TGeoPatternFinder *finder; //--- finder to be attached
1810 TString opt = ""; //--- option to be attached
1811 Double_t dphi;
1812 Int_t id;
1813 Double_t end = start+ndiv*step;
1814 switch (iaxis) {
1815 case 1: //--- R division
1816 finder = new TGeoPatternCylR(voldiv, ndiv, start, end);
1817 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1818 voldiv->SetFinder(finder);
1819 finder->SetDivIndex(voldiv->GetNdaughters());
1820 for (id=0; id<ndiv; id++) {
1821 shape = new TGeoTubeSeg(start+id*step, start+(id+1)*step, fDz, fPhi1, fPhi2);
1822 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1823 vmulti->AddVolume(vol);
1824 opt = "R";
1825 voldiv->AddNodeOffset(vol, id, 0, opt.Data());
1826 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1827 }
1828 return vmulti;
1829 case 2: //--- Phi division
1830 dphi = fPhi2-fPhi1;
1831 if (dphi<0) dphi+=360.;
1832 if (step<=0) {step=dphi/ndiv; start=fPhi1; end=fPhi2;}
1833 finder = new TGeoPatternCylPhi(voldiv, ndiv, start, end);
1834 voldiv->SetFinder(finder);
1835 finder->SetDivIndex(voldiv->GetNdaughters());
1836 shape = new TGeoTubeSeg(fRmin, fRmax, fDz, -step/2, step/2);
1837 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1838 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1839 vmulti->AddVolume(vol);
1840 opt = "Phi";
1841 for (id=0; id<ndiv; id++) {
1842 voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
1843 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1844 }
1845 return vmulti;
1846 case 3: //--- Z division
1847 finder = new TGeoPatternZ(voldiv, ndiv, start, end);
1848 voldiv->SetFinder(finder);
1849 finder->SetDivIndex(voldiv->GetNdaughters());
1850 shape = new TGeoTubeSeg(fRmin, fRmax, step/2, fPhi1, fPhi2);
1851 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1852 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1853 vmulti->AddVolume(vol);
1854 opt = "Z";
1855 for (id=0; id<ndiv; id++) {
1856 voldiv->AddNodeOffset(vol, id, start+step/2+id*step, opt.Data());
1857 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1858 }
1859 return vmulti;
1860 default:
1861 Error("Divide", "In shape %s wrong axis type for division", GetName());
1862 return 0;
1863 }
1864}
1865
1866////////////////////////////////////////////////////////////////////////////////
1867/// Get range of shape for a given axis.
1868
1870{
1871 xlo = 0;
1872 xhi = 0;
1873 Double_t dx = 0;
1874 switch (iaxis) {
1875 case 1:
1876 xlo = fRmin;
1877 xhi = fRmax;
1878 dx = xhi-xlo;
1879 return dx;
1880 case 2:
1881 xlo = fPhi1;
1882 xhi = fPhi2;
1883 dx = xhi-xlo;
1884 return dx;
1885 case 3:
1886 xlo = -fDz;
1887 xhi = fDz;
1888 dx = xhi-xlo;
1889 return dx;
1890 }
1891 return dx;
1892}
1893
1894////////////////////////////////////////////////////////////////////////////////
1895/// Fill vector param[4] with the bounding cylinder parameters. The order
1896/// is the following : Rmin, Rmax, Phi1, Phi2
1897
1899{
1900 param[0] = fRmin;
1901 param[0] *= param[0];
1902 param[1] = fRmax;
1903 param[1] *= param[1];
1904 param[2] = fPhi1;
1905 param[3] = fPhi2;
1906}
1907
1908////////////////////////////////////////////////////////////////////////////////
1909/// in case shape has some negative parameters, these has to be computed
1910/// in order to fit the mother
1911
1913{
1914 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
1915 if (!mother->TestShapeBit(kGeoTube)) {
1916 Error("GetMakeRuntimeShape", "Invalid mother for shape %s", GetName());
1917 return 0;
1918 }
1919 Double_t rmin, rmax, dz;
1920 rmin = fRmin;
1921 rmax = fRmax;
1922 dz = fDz;
1923 if (fDz<0) dz=((TGeoTube*)mother)->GetDz();
1924 if (fRmin<0)
1925 rmin = ((TGeoTube*)mother)->GetRmin();
1926 if ((fRmax<0) || (fRmax<=fRmin))
1927 rmax = ((TGeoTube*)mother)->GetRmax();
1928
1929 return (new TGeoTubeSeg(GetName(),rmin, rmax, dz, fPhi1, fPhi2));
1930}
1931
1932////////////////////////////////////////////////////////////////////////////////
1933/// print shape parameters
1934
1936{
1937 printf("*** Shape %s: TGeoTubeSeg ***\n", GetName());
1938 printf(" Rmin = %11.5f\n", fRmin);
1939 printf(" Rmax = %11.5f\n", fRmax);
1940 printf(" dz = %11.5f\n", fDz);
1941 printf(" phi1 = %11.5f\n", fPhi1);
1942 printf(" phi2 = %11.5f\n", fPhi2);
1943 printf(" Bounding box:\n");
1945}
1946
1947////////////////////////////////////////////////////////////////////////////////
1948/// Creates a TBuffer3D describing *this* shape.
1949/// Coordinates are in local reference frame.
1950
1952{
1954 Int_t nbPnts = 4*n;
1955 Int_t nbSegs = 2*nbPnts;
1956 Int_t nbPols = nbPnts-2;
1957
1959 nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols);
1960 if (buff)
1961 {
1962 SetPoints(buff->fPnts);
1963 SetSegsAndPols(*buff);
1964 }
1965
1966 return buff;
1967}
1968
1969////////////////////////////////////////////////////////////////////////////////
1970/// Fill TBuffer3D structure for segments and polygons.
1971
1973{
1974 Int_t i, j;
1976 Int_t c = GetBasicColor();
1977
1978 memset(buff.fSegs, 0, buff.NbSegs()*3*sizeof(Int_t));
1979 for (i = 0; i < 4; i++) {
1980 for (j = 1; j < n; j++) {
1981 buff.fSegs[(i*n+j-1)*3 ] = c;
1982 buff.fSegs[(i*n+j-1)*3+1] = i*n+j-1;
1983 buff.fSegs[(i*n+j-1)*3+2] = i*n+j;
1984 }
1985 }
1986 for (i = 4; i < 6; i++) {
1987 for (j = 0; j < n; j++) {
1988 buff.fSegs[(i*n+j)*3 ] = c+1;
1989 buff.fSegs[(i*n+j)*3+1] = (i-4)*n+j;
1990 buff.fSegs[(i*n+j)*3+2] = (i-2)*n+j;
1991 }
1992 }
1993 for (i = 6; i < 8; i++) {
1994 for (j = 0; j < n; j++) {
1995 buff.fSegs[(i*n+j)*3 ] = c;
1996 buff.fSegs[(i*n+j)*3+1] = 2*(i-6)*n+j;
1997 buff.fSegs[(i*n+j)*3+2] = (2*(i-6)+1)*n+j;
1998 }
1999 }
2000
2001 Int_t indx = 0;
2002 memset(buff.fPols, 0, buff.NbPols()*6*sizeof(Int_t));
2003 i = 0;
2004 for (j = 0; j < n-1; j++) {
2005 buff.fPols[indx++] = c;
2006 buff.fPols[indx++] = 4;
2007 buff.fPols[indx++] = (4+i)*n+j+1;
2008 buff.fPols[indx++] = (2+i)*n+j;
2009 buff.fPols[indx++] = (4+i)*n+j;
2010 buff.fPols[indx++] = i*n+j;
2011 }
2012 i = 1;
2013 for (j = 0; j < n-1; j++) {
2014 buff.fPols[indx++] = c;
2015 buff.fPols[indx++] = 4;
2016 buff.fPols[indx++] = i*n+j;
2017 buff.fPols[indx++] = (4+i)*n+j;
2018 buff.fPols[indx++] = (2+i)*n+j;
2019 buff.fPols[indx++] = (4+i)*n+j+1;
2020 }
2021 i = 2;
2022 for (j = 0; j < n-1; j++) {
2023 buff.fPols[indx++] = c+i;
2024 buff.fPols[indx++] = 4;
2025 buff.fPols[indx++] = (i-2)*2*n+j;
2026 buff.fPols[indx++] = (4+i)*n+j;
2027 buff.fPols[indx++] = ((i-2)*2+1)*n+j;
2028 buff.fPols[indx++] = (4+i)*n+j+1;
2029 }
2030 i = 3;
2031 for (j = 0; j < n-1; j++) {
2032 buff.fPols[indx++] = c+i;
2033 buff.fPols[indx++] = 4;
2034 buff.fPols[indx++] = (4+i)*n+j+1;
2035 buff.fPols[indx++] = ((i-2)*2+1)*n+j;
2036 buff.fPols[indx++] = (4+i)*n+j;
2037 buff.fPols[indx++] = (i-2)*2*n+j;
2038 }
2039 buff.fPols[indx++] = c+2;
2040 buff.fPols[indx++] = 4;
2041 buff.fPols[indx++] = 6*n;
2042 buff.fPols[indx++] = 4*n;
2043 buff.fPols[indx++] = 7*n;
2044 buff.fPols[indx++] = 5*n;
2045 buff.fPols[indx++] = c+2;
2046 buff.fPols[indx++] = 4;
2047 buff.fPols[indx++] = 6*n-1;
2048 buff.fPols[indx++] = 8*n-1;
2049 buff.fPols[indx++] = 5*n-1;
2050 buff.fPols[indx++] = 7*n-1;
2051}
2052
2053////////////////////////////////////////////////////////////////////////////////
2054/// computes the closest distance from given point InitTrigonometry();to this shape, according
2055/// to option. The matching point on the shape is stored in spoint.
2056
2058{
2059 Double_t safe = TGeoShape::Big();
2060 Double_t saf[3];
2061 Double_t rsq = point[0]*point[0]+point[1]*point[1];
2062 Double_t r = TMath::Sqrt(rsq);
2063 if (in) {
2064 saf[0] = fDz-TMath::Abs(point[2]);
2065 saf[1] = r-fRmin;
2066 saf[2] = fRmax-r;
2067 safe = saf[TMath::LocMin(3,saf)];
2068 if ((fPhi2-fPhi1)>=360.) return safe;
2069 Double_t safphi = TGeoShape::SafetyPhi(point,in,fPhi1,fPhi2);
2070 return TMath::Min(safe, safphi);
2071 }
2072 // Point expected to be outside
2073 Bool_t inphi = kFALSE;
2074 Double_t cpsi=point[0]*fCm+point[1]*fSm;
2075 saf[0] = TMath::Abs(point[2])-fDz;
2076 if (cpsi>r*fCdfi-TGeoShape::Tolerance()) inphi = kTRUE;
2077 if (inphi) {
2078 saf[1] = fRmin-r;
2079 saf[2] = r-fRmax;
2080 safe = saf[TMath::LocMax(3,saf)];
2081 safe = TMath::Max(0., safe);
2082 return safe;
2083 }
2084 // Point outside the phi range
2085 // Compute projected radius of the (r,phi) position vector onto
2086 // phi1 and phi2 edges and take the maximum for choosing the side.
2087 Double_t rproj = TMath::Max(point[0]*fC1+point[1]*fS1, point[0]*fC2+point[1]*fS2);
2088 saf[1] = fRmin-rproj;
2089 saf[2] = rproj-fRmax;
2090 safe = TMath::Max(saf[1], saf[2]);
2091 if ((fPhi2-fPhi1)>=360.) return TMath::Max(safe,saf[0]);
2092 if (safe>0) {
2093 // rproj not within (rmin,rmax) - > no need to calculate safphi
2094 safe = TMath::Sqrt(rsq-rproj*rproj+safe*safe);
2095 return (saf[0]<0) ? safe : TMath::Sqrt(safe*safe+saf[0]*saf[0]);
2096 }
2097 Double_t safphi = TGeoShape::SafetyPhi(point,in,fPhi1,fPhi2);
2098 return (saf[0]<0) ? safphi : TMath::Sqrt(saf[0]*saf[0]+safphi*safphi);
2099}
2100
2101////////////////////////////////////////////////////////////////////////////////
2102/// Static method to compute the closest distance from given point to this shape.
2103
2105 Double_t phi1d, Double_t phi2d, Int_t skipz)
2106{
2107 Double_t safe = TGeoShape::Big();
2108 Double_t saf[3];
2109 Double_t rsq = point[0]*point[0]+point[1]*point[1];
2110 Double_t r = TMath::Sqrt(rsq);
2111
2112 switch (skipz) {
2113 case 1: // skip lower Z plane
2114 saf[0] = dz - point[2];
2115 break;
2116 case 2: // skip upper Z plane
2117 saf[0] = dz + point[2];
2118 break;
2119 case 3: // skip both
2120 saf[0] = TGeoShape::Big();
2121 break;
2122 default:
2123 saf[0] = dz-TMath::Abs(point[2]);
2124 }
2125
2126 if (in) {
2127 saf[1] = r-rmin;
2128 saf[2] = rmax-r;
2129 safe = saf[TMath::LocMin(3,saf)];
2130 if ((phi2d-phi1d)>=360.) return safe;
2131 Double_t safphi = TGeoShape::SafetyPhi(point,in,phi1d,phi2d);
2132 return TMath::Min(safe, safphi);
2133 }
2134 // Point expected to be outside
2135 saf[0] = -saf[0];
2136 Bool_t inphi = kFALSE;
2137 Double_t phi1 = phi1d*TMath::DegToRad();
2138 Double_t phi2 = phi2d*TMath::DegToRad();
2139
2140 Double_t fio = 0.5*(phi1+phi2);
2141 Double_t cm = TMath::Cos(fio);
2142 Double_t sm = TMath::Sin(fio);
2143 Double_t cpsi=point[0]*cm+point[1]*sm;
2144 Double_t dfi = 0.5*(phi2-phi1);
2145 Double_t cdfi = TMath::Cos(dfi);
2146 if (cpsi>r*cdfi-TGeoShape::Tolerance()) inphi = kTRUE;
2147 if (inphi) {
2148 saf[1] = rmin-r;
2149 saf[2] = r-rmax;
2150 safe = saf[TMath::LocMax(3,saf)];
2151 safe = TMath::Max(0., safe);
2152 return safe;
2153 }
2154 // Point outside the phi range
2155 // Compute projected radius of the (r,phi) position vector onto
2156 // phi1 and phi2 edges and take the maximum for choosing the side.
2157 Double_t c1 = TMath::Cos(phi1);
2158 Double_t s1 = TMath::Sin(phi1);
2159 Double_t c2 = TMath::Cos(phi2);
2160 Double_t s2 = TMath::Sin(phi2);
2161
2162 Double_t rproj = TMath::Max(point[0]*c1+point[1]*s1, point[0]*c2+point[1]*s2);
2163 saf[1] = rmin-rproj;
2164 saf[2] = rproj-rmax;
2165 safe = TMath::Max(saf[1], saf[2]);
2166 if ((phi2d-phi1d)>=360.) return TMath::Max(safe,saf[0]);
2167 if (safe>0) {
2168 // rproj not within (rmin,rmax) - > no need to calculate safphi
2169 safe = TMath::Sqrt(rsq-rproj*rproj+safe*safe);
2170 return (saf[0]<0) ? safe : TMath::Sqrt(safe*safe+saf[0]*saf[0]);
2171 }
2172 Double_t safphi = TGeoShape::SafetyPhi(point,in,phi1d,phi2d);
2173 return (saf[0]<0) ? safphi : TMath::Sqrt(saf[0]*saf[0]+safphi*safphi);
2174}
2175
2176////////////////////////////////////////////////////////////////////////////////
2177/// Save a primitive as a C++ statement(s) on output stream "out".
2178
2179void TGeoTubeSeg::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
2180{
2182 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
2183 out << " rmin = " << fRmin << ";" << std::endl;
2184 out << " rmax = " << fRmax << ";" << std::endl;
2185 out << " dz = " << fDz << ";" << std::endl;
2186 out << " phi1 = " << fPhi1 << ";" << std::endl;
2187 out << " phi2 = " << fPhi2 << ";" << std::endl;
2188 out << " TGeoShape *" << GetPointerName() << " = new TGeoTubeSeg(\"" << GetName() << "\",rmin,rmax,dz,phi1,phi2);" << std::endl;
2190}
2191
2192////////////////////////////////////////////////////////////////////////////////
2193/// Set dimensions of the tube segment.
2194/// The segment will be from phiStart to phiEnd expressed in degree.
2195
2197 Double_t phiStart, Double_t phiEnd)
2198{
2199 fRmin = rmin;
2200 fRmax = rmax;
2201 fDz = dz;
2202 fPhi1 = phiStart;
2203 if (fPhi1 < 0) fPhi1 += 360.;
2204 fPhi2 = phiEnd;
2205 while (fPhi2<=fPhi1) fPhi2+=360.;
2206 if (TGeoShape::IsSameWithinTolerance(fPhi1,fPhi2)) Fatal("SetTubsDimensions", "In shape %s invalid phi1=%g, phi2=%g\n", GetName(), fPhi1, fPhi2);
2208}
2209
2210////////////////////////////////////////////////////////////////////////////////
2211/// Set dimensions of the tube segment starting from a list.
2212
2214{
2215 Double_t rmin = param[0];
2216 Double_t rmax = param[1];
2217 Double_t dz = param[2];
2218 Double_t phi1 = param[3];
2219 Double_t phi2 = param[4];
2220 SetTubsDimensions(rmin, rmax, dz, phi1, phi2);
2221}
2222
2223////////////////////////////////////////////////////////////////////////////////
2224/// Fills array with n random points located on the line segments of the shape mesh.
2225/// The output array must be provided with a length of minimum 3*npoints. Returns
2226/// true if operation is implemented.
2227
2229{
2230 if (npoints > (npoints/2)*2) {
2231 Error("GetPointsOnSegments","Npoints must be even number");
2232 return kFALSE;
2233 }
2234 Int_t nc = (Int_t)TMath::Sqrt(0.5*npoints);
2235 Double_t dphi = (fPhi2-fPhi1)*TMath::DegToRad()/(nc-1);
2236 Double_t phi = 0;
2237 Double_t phi1 = fPhi1 * TMath::DegToRad();
2238 Int_t ntop = npoints/2 - nc*(nc-1);
2239 Double_t dz = 2*fDz/(nc-1);
2240 Double_t z = 0;
2241 Int_t icrt = 0;
2242 Int_t nphi = nc;
2243 // loop z sections
2244 for (Int_t i=0; i<nc; i++) {
2245 if (i == (nc-1)) {
2246 nphi = ntop;
2247 dphi = (fPhi2-fPhi1)*TMath::DegToRad()/(nphi-1);
2248 }
2249 z = -fDz + i*dz;
2250 // loop points on circle sections
2251 for (Int_t j=0; j<nphi; j++) {
2252 phi = phi1 + j*dphi;
2253 array[icrt++] = fRmin * TMath::Cos(phi);
2254 array[icrt++] = fRmin * TMath::Sin(phi);
2255 array[icrt++] = z;
2256 array[icrt++] = fRmax * TMath::Cos(phi);
2257 array[icrt++] = fRmax * TMath::Sin(phi);
2258 array[icrt++] = z;
2259 }
2260 }
2261 return kTRUE;
2262}
2263
2264////////////////////////////////////////////////////////////////////////////////
2265/// Create tube segment mesh points.
2266
2268{
2269 Double_t dz;
2270 Int_t j, n;
2271 Double_t phi, phi1, phi2, dphi;
2272 phi1 = fPhi1;
2273 phi2 = fPhi2;
2274 if (phi2<phi1) phi2+=360.;
2275 n = gGeoManager->GetNsegments()+1;
2276
2277 dphi = (phi2-phi1)/(n-1);
2278 dz = fDz;
2279
2280 if (points) {
2281 Int_t indx = 0;
2282
2283 for (j = 0; j < n; j++) {
2284 phi = (phi1+j*dphi)*TMath::DegToRad();
2285 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
2286 indx++;
2287 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
2288 indx++;
2289 points[indx+6*n] = dz;
2290 points[indx] =-dz;
2291 indx++;
2292 }
2293 for (j = 0; j < n; j++) {
2294 phi = (phi1+j*dphi)*TMath::DegToRad();
2295 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
2296 indx++;
2297 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
2298 indx++;
2299 points[indx+6*n]= dz;
2300 points[indx] =-dz;
2301 indx++;
2302 }
2303 }
2304}
2305
2306////////////////////////////////////////////////////////////////////////////////
2307/// Create tube segment mesh points.
2308
2310{
2311 Double_t dz;
2312 Int_t j, n;
2313 Double_t phi, phi1, phi2, dphi;
2314 phi1 = fPhi1;
2315 phi2 = fPhi2;
2316 if (phi2<phi1) phi2+=360.;
2317 n = gGeoManager->GetNsegments()+1;
2318
2319 dphi = (phi2-phi1)/(n-1);
2320 dz = fDz;
2321
2322 if (points) {
2323 Int_t indx = 0;
2324
2325 for (j = 0; j < n; j++) {
2326 phi = (phi1+j*dphi)*TMath::DegToRad();
2327 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
2328 indx++;
2329 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
2330 indx++;
2331 points[indx+6*n] = dz;
2332 points[indx] =-dz;
2333 indx++;
2334 }
2335 for (j = 0; j < n; j++) {
2336 phi = (phi1+j*dphi)*TMath::DegToRad();
2337 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
2338 indx++;
2339 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
2340 indx++;
2341 points[indx+6*n]= dz;
2342 points[indx] =-dz;
2343 indx++;
2344 }
2345 }
2346}
2347
2348////////////////////////////////////////////////////////////////////////////////
2349/// Returns numbers of vertices, segments and polygons composing the shape mesh.
2350
2351void TGeoTubeSeg::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
2352{
2354 nvert = n*4;
2355 nsegs = n*8;
2356 npols = n*4 - 2;
2357}
2358
2359////////////////////////////////////////////////////////////////////////////////
2360/// Return number of vertices of the mesh representation
2361
2363{
2365 Int_t numPoints = n*4;
2366 return numPoints;
2367}
2368
2369////////////////////////////////////////////////////////////////////////////////
2370/// fill size of this 3-D object
2371
2373{
2374}
2375
2376////////////////////////////////////////////////////////////////////////////////
2377/// Fills a static 3D buffer and returns a reference.
2378
2379const TBuffer3D & TGeoTubeSeg::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
2380{
2381 static TBuffer3DTubeSeg buffer;
2382 TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
2383
2384 if (reqSections & TBuffer3D::kShapeSpecific) {
2385 // These from TBuffer3DTube / TGeoTube
2386 buffer.fRadiusInner = fRmin;
2387 buffer.fRadiusOuter = fRmax;
2388 buffer.fHalfLength = fDz;
2389 buffer.fPhiMin = fPhi1;
2390 buffer.fPhiMax = fPhi2;
2392 }
2393 if (reqSections & TBuffer3D::kRawSizes) {
2395 Int_t nbPnts = 4*n;
2396 Int_t nbSegs = 2*nbPnts;
2397 Int_t nbPols = nbPnts-2;
2398 if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
2400 }
2401 }
2402 if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
2403 SetPoints(buffer.fPnts);
2404 if (!buffer.fLocalFrame) {
2405 TransformPoints(buffer.fPnts, buffer.NbPnts());
2406 }
2407 SetSegsAndPols(buffer);
2409 }
2410
2411 return buffer;
2412}
2413
2414////////////////////////////////////////////////////////////////////////////////
2415/// Check the inside status for each of the points in the array.
2416/// Input: Array of point coordinates + vector size
2417/// Output: Array of Booleans for the inside of each point
2418
2419void TGeoTubeSeg::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
2420{
2421 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
2422}
2423
2424////////////////////////////////////////////////////////////////////////////////
2425/// Compute the normal for an array o points so that norm.dot.dir is positive
2426/// Input: Arrays of point coordinates and directions + vector size
2427/// Output: Array of normal directions
2428
2429void TGeoTubeSeg::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
2430{
2431 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
2432}
2433
2434////////////////////////////////////////////////////////////////////////////////
2435/// Compute distance from array of input points having directions specified by dirs. Store output in dists
2436
2437void TGeoTubeSeg::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2438{
2439 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
2440}
2441
2442////////////////////////////////////////////////////////////////////////////////
2443/// Compute distance from array of input points having directions specified by dirs. Store output in dists
2444
2445void TGeoTubeSeg::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2446{
2447 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
2448}
2449
2450////////////////////////////////////////////////////////////////////////////////
2451/// Compute safe distance from each of the points in the input array.
2452/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
2453/// Output: Safety values
2454
2455void TGeoTubeSeg::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
2456{
2457 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
2458}
2459
2461
2463{
2464// default ctor
2465 fNlow[0] = fNlow[1] = fNhigh[0] = fNhigh[1] = 0.;
2466 fNlow[2] = -1;
2467 fNhigh[2] = 1;
2468}
2469
2470////////////////////////////////////////////////////////////////////////////////
2471/// constructor
2472
2474 Double_t lx, Double_t ly, Double_t lz, Double_t tx, Double_t ty, Double_t tz)
2475 :TGeoTubeSeg(rmin, rmax, dz, phi1, phi2)
2476{
2477 fNlow[0] = lx;
2478 fNlow[1] = ly;
2479 fNlow[2] = lz;
2480 fNhigh[0] = tx;
2481 fNhigh[1] = ty;
2482 fNhigh[2] = tz;
2484 ComputeBBox();
2485}
2486
2487////////////////////////////////////////////////////////////////////////////////
2488/// constructor
2489
2490TGeoCtub::TGeoCtub(const char *name, Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2,
2491 Double_t lx, Double_t ly, Double_t lz, Double_t tx, Double_t ty, Double_t tz)
2492 :TGeoTubeSeg(name, rmin, rmax, dz, phi1, phi2)
2493{
2494 fNlow[0] = lx;
2495 fNlow[1] = ly;
2496 fNlow[2] = lz;
2497 fNhigh[0] = tx;
2498 fNhigh[1] = ty;
2499 fNhigh[2] = tz;
2501 ComputeBBox();
2502}
2503
2504////////////////////////////////////////////////////////////////////////////////
2505/// ctor with parameters
2506
2508 :TGeoTubeSeg(0,0,0,0,0)
2509{
2510 SetCtubDimensions(params[0], params[1], params[2], params[3], params[4], params[5],
2511 params[6], params[7], params[8], params[9], params[10]);
2513}
2514
2515////////////////////////////////////////////////////////////////////////////////
2516/// destructor
2517
2519{
2520}
2521
2522////////////////////////////////////////////////////////////////////////////////
2523/// Computes capacity of the shape in [length^3]
2524
2526{
2527 Double_t capacity = TGeoTubeSeg::Capacity();
2528 return capacity;
2529}
2530
2531////////////////////////////////////////////////////////////////////////////////
2532/// compute minimum bounding box of the ctub
2533
2535{
2537 if ((fNlow[2]>-(1E-10)) || (fNhigh[2]<1E-10)) {
2538 Error("ComputeBBox", "In shape %s wrong definition of cut planes", GetName());
2539 return;
2540 }
2541 Double_t xc=0, yc=0;
2542 Double_t zmin=0, zmax=0;
2543 Double_t z1;
2544 Double_t z[8];
2545 // check if nxy is in the phi range
2546 Double_t phi_low = TMath::ATan2(fNlow[1], fNlow[0]) *TMath::RadToDeg();
2548 Bool_t in_range_low = kFALSE;
2549 Bool_t in_range_hi = kFALSE;
2550
2551 Int_t i;
2552 for (i=0; i<2; i++) {
2553 if (phi_low<0) phi_low+=360.;
2554 Double_t dphi = fPhi2 -fPhi1;
2555 if (dphi < 0) dphi+=360.;
2556 Double_t ddp = phi_low-fPhi1;
2557 if (ddp<0) ddp += 360.;
2558 if (ddp <= dphi) {
2559 xc = fRmin*TMath::Cos(phi_low*TMath::DegToRad());
2560 yc = fRmin*TMath::Sin(phi_low*TMath::DegToRad());
2561 z1 = GetZcoord(xc, yc, -fDz);
2562 xc = fRmax*TMath::Cos(phi_low*TMath::DegToRad());
2563 yc = fRmax*TMath::Sin(phi_low*TMath::DegToRad());
2564 z1 = TMath::Min(z1, GetZcoord(xc, yc, -fDz));
2565 if (in_range_low)
2566 zmin = TMath::Min(zmin, z1);
2567 else
2568 zmin = z1;
2569 in_range_low = kTRUE;
2570 }
2571 phi_low += 180;
2572 if (phi_low>360) phi_low-=360.;
2573 }
2574
2575 for (i=0; i<2; i++) {
2576 if (phi_hi<0) phi_hi+=360.;
2577 Double_t dphi = fPhi2 -fPhi1;
2578 if (dphi < 0) dphi+=360.;
2579 Double_t ddp = phi_hi-fPhi1;
2580 if (ddp<0) ddp += 360.;
2581 if (ddp <= dphi) {
2582 xc = fRmin*TMath::Cos(phi_hi*TMath::DegToRad());
2583 yc = fRmin*TMath::Sin(phi_hi*TMath::DegToRad());
2584 z1 = GetZcoord(xc, yc, fDz);
2585 xc = fRmax*TMath::Cos(phi_hi*TMath::DegToRad());
2586 yc = fRmax*TMath::Sin(phi_hi*TMath::DegToRad());
2587 z1 = TMath::Max(z1, GetZcoord(xc, yc, fDz));
2588 if (in_range_hi)
2589 zmax = TMath::Max(zmax, z1);
2590 else
2591 zmax = z1;
2592 in_range_hi = kTRUE;
2593 }
2594 phi_hi += 180;
2595 if (phi_hi>360) phi_hi-=360.;
2596 }
2597
2598
2599 xc = fRmin*fC1;
2600 yc = fRmin*fS1;
2601 z[0] = GetZcoord(xc, yc, -fDz);
2602 z[4] = GetZcoord(xc, yc, fDz);
2603
2604 xc = fRmin*fC2;
2605 yc = fRmin*fS2;
2606 z[1] = GetZcoord(xc, yc, -fDz);
2607 z[5] = GetZcoord(xc, yc, fDz);
2608
2609 xc = fRmax*fC1;
2610 yc = fRmax*fS1;
2611 z[2] = GetZcoord(xc, yc, -fDz);
2612 z[6] = GetZcoord(xc, yc, fDz);
2613
2614 xc = fRmax*fC2;
2615 yc = fRmax*fS2;
2616 z[3] = GetZcoord(xc, yc, -fDz);
2617 z[7] = GetZcoord(xc, yc, fDz);
2618
2619 z1 = z[TMath::LocMin(4, &z[0])];
2620 if (in_range_low)
2621 zmin = TMath::Min(zmin, z1);
2622 else
2623 zmin = z1;
2624
2625 z1 = z[TMath::LocMax(4, &z[4])+4];
2626 if (in_range_hi)
2627 zmax = TMath::Max(zmax, z1);
2628 else
2629 zmax = z1;
2630
2631 fDZ = 0.5*(zmax-zmin);
2632 fOrigin[2] = 0.5*(zmax+zmin);
2633}
2634
2635////////////////////////////////////////////////////////////////////////////////
2636/// Compute normal to closest surface from POINT.
2637
2638void TGeoCtub::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
2639{
2640 Double_t saf[4];
2641 Bool_t isseg = kTRUE;
2642 if (TMath::Abs(fPhi2-fPhi1-360.)<1E-8) isseg=kFALSE;
2643 Double_t rsq = point[0]*point[0]+point[1]*point[1];
2644 Double_t r = TMath::Sqrt(rsq);
2645
2646 saf[0] = TMath::Abs(point[0]*fNlow[0] + point[1]*fNlow[1] + (fDz+point[2])*fNlow[2]);
2647 saf[1] = TMath::Abs(point[0]*fNhigh[0] + point[1]*fNhigh[1] - (fDz-point[2])*fNhigh[2]);
2648 saf[2] = (fRmin>1E-10)?TMath::Abs(r-fRmin):TGeoShape::Big();
2649 saf[3] = TMath::Abs(fRmax-r);
2650 Int_t i = TMath::LocMin(4,saf);
2651 if (isseg) {
2652 if (TGeoShape::IsCloseToPhi(saf[i], point,fC1,fS1,fC2,fS2)) {
2653 TGeoShape::NormalPhi(point,dir,norm,fC1,fS1,fC2,fS2);
2654 return;
2655 }
2656 }
2657 if (i==0) {
2658 memcpy(norm, fNlow, 3*sizeof(Double_t));
2659 if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
2660 norm[0] = -norm[0];
2661 norm[1] = -norm[1];
2662 norm[2] = -norm[2];
2663 }
2664 return;
2665 }
2666 if (i==1) {
2667 memcpy(norm, fNhigh, 3*sizeof(Double_t));
2668 if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
2669 norm[0] = -norm[0];
2670 norm[1] = -norm[1];
2671 norm[2] = -norm[2];
2672 }
2673 return;
2674 }
2675
2676 norm[2] = 0;
2677 Double_t phi = TMath::ATan2(point[1], point[0]);
2678 norm[0] = TMath::Cos(phi);
2679 norm[1] = TMath::Sin(phi);
2680 if (norm[0]*dir[0]+norm[1]*dir[1]<0) {
2681 norm[0] = -norm[0];
2682 norm[1] = -norm[1];
2683 }
2684}
2685
2686////////////////////////////////////////////////////////////////////////////////
2687/// check if point is contained in the cut tube
2688/// check the lower cut plane
2689
2691{
2692 Double_t zin = point[0]*fNlow[0]+point[1]*fNlow[1]+(point[2]+fDz)*fNlow[2];
2693 if (zin>0) return kFALSE;
2694 // check the higher cut plane
2695 zin = point[0]*fNhigh[0]+point[1]*fNhigh[1]+(point[2]-fDz)*fNhigh[2];
2696 if (zin>0) return kFALSE;
2697 // check radius
2698 Double_t r2 = point[0]*point[0]+point[1]*point[1];
2699 if ((r2<fRmin*fRmin) || (r2>fRmax*fRmax)) return kFALSE;
2700 // check phi
2701 Double_t phi = TMath::ATan2(point[1], point[0]) * TMath::RadToDeg();
2702 if (phi < 0 ) phi+=360.;
2703 Double_t dphi = fPhi2 -fPhi1;
2704 Double_t ddp = phi-fPhi1;
2705 if (ddp<0) ddp += 360.;
2706// if (ddp>360) ddp-=360;
2707 if (ddp > dphi) return kFALSE;
2708 return kTRUE;
2709}
2710
2711////////////////////////////////////////////////////////////////////////////////
2712/// Get range of shape for a given axis.
2713
2715{
2716 xlo = 0;
2717 xhi = 0;
2718 Double_t dx = 0;
2719 switch (iaxis) {
2720 case 1:
2721 xlo = fRmin;
2722 xhi = fRmax;
2723 dx = xhi-xlo;
2724 return dx;
2725 case 2:
2726 xlo = fPhi1;
2727 xhi = fPhi2;
2728 dx = xhi-xlo;
2729 return dx;
2730 }
2731 return dx;
2732}
2733
2734////////////////////////////////////////////////////////////////////////////////
2735/// compute real Z coordinate of a point belonging to either lower or
2736/// higher caps (z should be either +fDz or -fDz)
2737
2739{
2740 Double_t newz = 0;
2741 if (zc<0) newz = -fDz-(xc*fNlow[0]+yc*fNlow[1])/fNlow[2];
2742 else newz = fDz-(xc*fNhigh[0]+yc*fNhigh[1])/fNhigh[2];
2743 return newz;
2744}
2745
2746////////////////////////////////////////////////////////////////////////////////
2747/// compute distance from outside point to surface of the cut tube
2748
2749Double_t TGeoCtub::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
2750{
2751 if (iact<3 && safe) {
2752 *safe = Safety(point, kFALSE);
2753 if (iact==0) return TGeoShape::Big();
2754 if ((iact==1) && (step<=*safe)) return TGeoShape::Big();
2755 }
2756// Check if the bounding box is crossed within the requested distance
2757 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
2758 if (sdist>=step) return TGeoShape::Big();
2759 Double_t saf[2];
2760 saf[0] = point[0]*fNlow[0] + point[1]*fNlow[1] + (fDz+point[2])*fNlow[2];
2761 saf[1] = point[0]*fNhigh[0] + point[1]*fNhigh[1] + (point[2]-fDz)*fNhigh[2];
2762 Double_t rsq = point[0]*point[0]+point[1]*point[1];
2763 Double_t r = TMath::Sqrt(rsq);
2764 Double_t cpsi=0;
2765 Bool_t tub = kFALSE;
2766 if (TMath::Abs(fPhi2-fPhi1-360.)<1E-8) tub = kTRUE;
2767
2768 // find distance to shape
2769 Double_t r2;
2770 Double_t calf = dir[0]*fNlow[0]+dir[1]*fNlow[1]+dir[2]*fNlow[2];
2771 // check Z planes
2772 Double_t xi, yi, zi;
2774 if (saf[0]>0) {
2775 if (calf<0) {
2776 s = -saf[0]/calf;
2777 xi = point[0]+s*dir[0];
2778 yi = point[1]+s*dir[1];
2779 r2=xi*xi+yi*yi;
2780 if (((fRmin*fRmin)<=r2) && (r2<=(fRmax*fRmax))) {
2781 if (tub) return s;
2782 cpsi=(xi*fCm+yi*fSm)/TMath::Sqrt(r2);
2783 if (cpsi>=fCdfi) return s;
2784 }
2785 }
2786 }
2787 calf = dir[0]*fNhigh[0]+dir[1]*fNhigh[1]+dir[2]*fNhigh[2];
2788 if (saf[1]>0) {
2789 if (calf<0) {
2790 s = -saf[1]/calf;
2791 xi = point[0]+s*dir[0];
2792 yi = point[1]+s*dir[1];
2793 r2=xi*xi+yi*yi;
2794 if (((fRmin*fRmin)<=r2) && (r2<=(fRmax*fRmax))) {
2795 if (tub) return s;
2796 cpsi=(xi*fCm+yi*fSm)/TMath::Sqrt(r2);
2797 if (cpsi>=fCdfi) return s;
2798 }
2799 }
2800 }
2801
2802 // check outer cyl. surface
2803 Double_t nsq=dir[0]*dir[0]+dir[1]*dir[1];
2804 if (TMath::Abs(nsq)<1E-10) return TGeoShape::Big();
2805 Double_t rdotn=point[0]*dir[0]+point[1]*dir[1];
2806 Double_t b,d;
2807 // only r>fRmax coming inwards has to be considered
2808 if (r>fRmax && rdotn<0) {
2809 TGeoTube::DistToTube(rsq, nsq, rdotn, fRmax, b, d);
2810 if (d>0) {
2811 s=-b-d;
2812 if (s>0) {
2813 xi=point[0]+s*dir[0];
2814 yi=point[1]+s*dir[1];
2815 zi=point[2]+s*dir[2];
2816 if ((-xi*fNlow[0]-yi*fNlow[1]-(zi+fDz)*fNlow[2])>0) {
2817 if ((-xi*fNhigh[0]-yi*fNhigh[1]+(fDz-zi)*fNhigh[2])>0) {
2818 if (tub) return s;
2819 cpsi=(xi*fCm+yi*fSm)/fRmax;
2820 if (cpsi>=fCdfi) return s;
2821 }
2822 }
2823 }
2824 }
2825 }
2826 // check inner cylinder
2827 Double_t snxt=TGeoShape::Big();
2828 if (fRmin>0) {
2829 TGeoTube::DistToTube(rsq, nsq, rdotn, fRmin, b, d);
2830 if (d>0) {
2831 s=-b+d;
2832 if (s>0) {
2833 xi=point[0]+s*dir[0];
2834 yi=point[1]+s*dir[1];
2835 zi=point[2]+s*dir[2];
2836 if ((-xi*fNlow[0]-yi*fNlow[1]-(zi+fDz)*fNlow[2])>0) {
2837 if ((-xi*fNhigh[0]-yi*fNhigh[1]+(fDz-zi)*fNhigh[2])>0) {
2838 if (tub) return s;
2839 cpsi=(xi*fCm+yi*fSm)/fRmin;
2840 if (cpsi>=fCdfi) snxt=s;
2841 }
2842 }
2843 }
2844 }
2845 }
2846 // check phi planes
2847 if (tub) return snxt;
2848 Double_t un=dir[0]*fS1-dir[1]*fC1;
2849 if (un<-TGeoShape::Tolerance()) {
2850 s=(point[1]*fC1-point[0]*fS1)/un;
2851 if (s>=0) {
2852 xi=point[0]+s*dir[0];
2853 yi=point[1]+s*dir[1];
2854 zi=point[2]+s*dir[2];
2855 if ((-xi*fNlow[0]-yi*fNlow[1]-(zi+fDz)*fNlow[2])>0) {
2856 if ((-xi*fNhigh[0]-yi*fNhigh[1]+(fDz-zi)*fNhigh[2])>0) {
2857 r2=xi*xi+yi*yi;
2858 if ((fRmin*fRmin<=r2) && (r2<=fRmax*fRmax)) {
2859 if ((yi*fCm-xi*fSm)<=0) {
2860 if (s<snxt) snxt=s;
2861 }
2862 }
2863 }
2864 }
2865 }
2866 }
2867 un=dir[0]*fS2-dir[1]*fC2;
2868 if (un>TGeoShape::Tolerance()) {
2869 s=(point[1]*fC2-point[0]*fS2)/un;
2870 if (s>=0) {
2871 xi=point[0]+s*dir[0];
2872 yi=point[1]+s*dir[1];
2873 zi=point[2]+s*dir[2];
2874 if ((-xi*fNlow[0]-yi*fNlow[1]-(zi+fDz)*fNlow[2])>0) {
2875 if ((-xi*fNhigh[0]-yi*fNhigh[1]+(fDz-zi)*fNhigh[2])>0) {
2876 r2=xi*xi+yi*yi;
2877 if ((fRmin*fRmin<=r2) && (r2<=fRmax*fRmax)) {
2878 if ((yi*fCm-xi*fSm)>=0) {
2879 if (s<snxt) snxt=s;
2880 }
2881 }
2882 }
2883 }
2884 }
2885 }
2886 return snxt;
2887}
2888
2889////////////////////////////////////////////////////////////////////////////////
2890/// compute distance from inside point to surface of the cut tube
2891
2892Double_t TGeoCtub::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
2893{
2894 if (iact<3 && safe) *safe = Safety(point, kTRUE);
2895 if (iact==0) return TGeoShape::Big();
2896 if ((iact==1) && (*safe>step)) return TGeoShape::Big();
2897 Double_t rsq = point[0]*point[0]+point[1]*point[1];
2898 Bool_t tub = kFALSE;
2899 if (TMath::Abs(fPhi2-fPhi1-360.)<1E-8) tub = kTRUE;
2900 // compute distance to surface
2901 // Do Z
2902 Double_t sz = TGeoShape::Big();
2903 Double_t saf[2];
2904 saf[0] = -point[0]*fNlow[0] - point[1]*fNlow[1] - (fDz+point[2])*fNlow[2];
2905 saf[1] = -point[0]*fNhigh[0] - point[1]*fNhigh[1] + (fDz-point[2])*fNhigh[2];
2906 Double_t calf = dir[0]*fNlow[0]+dir[1]*fNlow[1]+dir[2]*fNlow[2];
2907 if (calf>0) sz = saf[0]/calf;
2908
2910 calf = dir[0]*fNhigh[0]+dir[1]*fNhigh[1]+dir[2]*fNhigh[2];
2911 if (calf>0) {
2912 sz1 = saf[1]/calf;
2913 if (sz1<sz) sz = sz1;
2914 }
2915
2916 // Do R
2917 Double_t nsq=dir[0]*dir[0]+dir[1]*dir[1];
2918 // track parallel to Z
2919 if (TMath::Abs(nsq)<1E-10) return sz;
2920 Double_t rdotn=point[0]*dir[0]+point[1]*dir[1];
2922 Double_t b, d;
2923 Bool_t skip_outer = kFALSE;
2924 // inner cylinder
2925 if (fRmin>1E-10) {
2926 TGeoTube::DistToTube(rsq, nsq, rdotn, fRmin, b, d);
2927 if (d>0) {
2928 sr=-b-d;
2929 if (sr>0) skip_outer = kTRUE;
2930 }
2931 }
2932 // outer cylinder
2933 if (!skip_outer) {
2934 TGeoTube::DistToTube(rsq, nsq, rdotn, fRmax, b, d);
2935 if (d>0) {
2936 sr=-b+d;
2937 if (sr<0) sr=TGeoShape::Big();
2938 } else {
2939 return 0.; // already outside
2940 }
2941 }
2942 // phi planes
2943 Double_t sfmin = TGeoShape::Big();
2944 if (!tub) sfmin=TGeoShape::DistToPhiMin(point, dir, fS1, fC1, fS2, fC2, fSm, fCm);
2945 return TMath::Min(TMath::Min(sz,sr), sfmin);
2946}
2947
2948////////////////////////////////////////////////////////////////////////////////
2949/// Divide the tube along one axis.
2950
2951TGeoVolume *TGeoCtub::Divide(TGeoVolume * /*voldiv*/, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/,
2952 Double_t /*start*/, Double_t /*step*/)
2953{
2954 Warning("Divide", "In shape %s division of a cut tube not implemented", GetName());
2955 return 0;
2956}
2957
2958////////////////////////////////////////////////////////////////////////////////
2959/// in case shape has some negative parameters, these has to be computed
2960/// in order to fit the mother
2961
2963{
2964 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
2965 if (!mother->TestShapeBit(kGeoTube)) {
2966 Error("GetMakeRuntimeShape", "Invalid mother for shape %s", GetName());
2967 return 0;
2968 }
2969 Double_t rmin, rmax, dz;
2970 rmin = fRmin;
2971 rmax = fRmax;
2972 dz = fDz;
2973 if (fDz<0) dz=((TGeoTube*)mother)->GetDz();
2974 if (fRmin<0)
2975 rmin = ((TGeoTube*)mother)->GetRmin();
2976 if ((fRmax<0) || (fRmax<=fRmin))
2977 rmax = ((TGeoTube*)mother)->GetRmax();
2978
2979 return (new TGeoCtub(rmin, rmax, dz, fPhi1, fPhi2, fNlow[0], fNlow[1], fNlow[2],
2980 fNhigh[0], fNhigh[1], fNhigh[2]));
2981}
2982
2983////////////////////////////////////////////////////////////////////////////////
2984/// print shape parameters
2985
2987{
2988 printf("*** Shape %s: TGeoCtub ***\n", GetName());
2989 printf(" lx = %11.5f\n", fNlow[0]);
2990 printf(" ly = %11.5f\n", fNlow[1]);
2991 printf(" lz = %11.5f\n", fNlow[2]);
2992 printf(" tx = %11.5f\n", fNhigh[0]);
2993 printf(" ty = %11.5f\n", fNhigh[1]);
2994 printf(" tz = %11.5f\n", fNhigh[2]);
2996}
2997
2998////////////////////////////////////////////////////////////////////////////////
2999/// computes the closest distance from given point to this shape, according
3000/// to option. The matching point on the shape is stored in spoint.
3001
3003{
3004 Double_t saf[4];
3005 Double_t rsq = point[0]*point[0]+point[1]*point[1];
3006 Double_t r = TMath::Sqrt(rsq);
3007 Bool_t isseg = kTRUE;
3008 if (TMath::Abs(fPhi2-fPhi1-360.)<1E-8) isseg=kFALSE;
3009
3010 saf[0] = -point[0]*fNlow[0] - point[1]*fNlow[1] - (fDz+point[2])*fNlow[2];
3011 saf[1] = -point[0]*fNhigh[0] - point[1]*fNhigh[1] + (fDz-point[2])*fNhigh[2];
3012 saf[2] = (fRmin<1E-10 && !isseg)?TGeoShape::Big():(r-fRmin);
3013 saf[3] = fRmax-r;
3014 Double_t safphi = TGeoShape::Big();
3015 Double_t safe = TGeoShape::Big();
3016 if (isseg) safphi = TGeoShape::SafetyPhi(point, in, fPhi1, fPhi2);
3017
3018 if (in) {
3019 safe = saf[TMath::LocMin(4,saf)];
3020 return TMath::Min(safe, safphi);
3021 }
3022 for (Int_t i=0; i<4; i++) saf[i]=-saf[i];
3023 safe = saf[TMath::LocMax(4,saf)];
3024 if (isseg) return TMath::Max(safe, safphi);
3025 return safe;
3026}
3027
3028////////////////////////////////////////////////////////////////////////////////
3029/// set dimensions of a cut tube
3030
3032 Double_t lx, Double_t ly, Double_t lz, Double_t tx, Double_t ty, Double_t tz)
3033{
3034 SetTubsDimensions(rmin, rmax, dz, phi1, phi2);
3035 fNlow[0] = lx;
3036 fNlow[1] = ly;
3037 fNlow[2] = lz;
3038 fNhigh[0] = tx;
3039 fNhigh[1] = ty;
3040 fNhigh[2] = tz;
3041 ComputeBBox();
3042}
3043
3044////////////////////////////////////////////////////////////////////////////////
3045/// Save a primitive as a C++ statement(s) on output stream "out".
3046
3047void TGeoCtub::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
3048{
3050 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
3051 out << " rmin = " << fRmin << ";" << std::endl;
3052 out << " rmax = " << fRmax << ";" << std::endl;
3053 out << " dz = " << fDz << ";" << std::endl;
3054 out << " phi1 = " << fPhi1 << ";" << std::endl;
3055 out << " phi2 = " << fPhi2 << ";" << std::endl;
3056 out << " lx = " << fNlow[0] << ";" << std::endl;
3057 out << " ly = " << fNlow[1] << ";" << std::endl;
3058 out << " lz = " << fNlow[2] << ";" << std::endl;
3059 out << " tx = " << fNhigh[0] << ";" << std::endl;
3060 out << " ty = " << fNhigh[1] << ";" << std::endl;
3061 out << " tz = " << fNhigh[2] << ";" << std::endl;
3062 out << " TGeoShape *" << GetPointerName() << " = new TGeoCtub(\"" << GetName() << "\",rmin,rmax,dz,phi1,phi2,lx,ly,lz,tx,ty,tz);" << std::endl; TObject::SetBit(TGeoShape::kGeoSavePrimitive);
3063}
3064
3065////////////////////////////////////////////////////////////////////////////////
3066/// Set dimensions of the cut tube starting from a list.
3067
3069{
3070 SetCtubDimensions(param[0], param[1], param[2], param[3], param[4], param[5],
3071 param[6], param[7], param[8], param[9], param[10]);
3072 ComputeBBox();
3073}
3074
3075////////////////////////////////////////////////////////////////////////////////
3076/// Fills array with n random points located on the line segments of the shape mesh.
3077/// The output array must be provided with a length of minimum 3*npoints. Returns
3078/// true if operation is implemented.
3079
3081{
3082 return kFALSE;
3083}
3084
3085////////////////////////////////////////////////////////////////////////////////
3086/// Create mesh points for the cut tube.
3087
3089{
3090 Double_t dz;
3091 Int_t j, n;
3092 Double_t phi, phi1, phi2, dphi;
3093 phi1 = fPhi1;
3094 phi2 = fPhi2;
3095 if (phi2<phi1) phi2+=360.;
3096 n = gGeoManager->GetNsegments()+1;
3097
3098 dphi = (phi2-phi1)/(n-1);
3099 dz = fDz;
3100
3101 if (points) {
3102 Int_t indx = 0;
3103
3104 for (j = 0; j < n; j++) {
3105 phi = (phi1+j*dphi)*TMath::DegToRad();
3106 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
3107 indx++;
3108 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
3109 indx++;
3110 points[indx+6*n] = GetZcoord(points[indx-2], points[indx-1], dz);
3111 points[indx] = GetZcoord(points[indx-2], points[indx-1], -dz);
3112 indx++;
3113 }
3114 for (j = 0; j < n; j++) {
3115 phi = (phi1+j*dphi)*TMath::DegToRad();
3116 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
3117 indx++;
3118 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
3119 indx++;
3120 points[indx+6*n]= GetZcoord(points[indx-2], points[indx-1], dz);
3121 points[indx] = GetZcoord(points[indx-2], points[indx-1], -dz);
3122 indx++;
3123 }
3124 }
3125}
3126
3127////////////////////////////////////////////////////////////////////////////////
3128/// Create mesh points for the cut tube.
3129
3131{
3132 Double_t dz;
3133 Int_t j, n;
3134 Double_t phi, phi1, phi2, dphi;
3135 phi1 = fPhi1;
3136 phi2 = fPhi2;
3137 if (phi2<phi1) phi2+=360.;
3138 n = gGeoManager->GetNsegments()+1;
3139
3140 dphi = (phi2-phi1)/(n-1);
3141 dz = fDz;
3142
3143 if (points) {
3144 Int_t indx = 0;
3145
3146 for (j = 0; j < n; j++) {
3147 phi = (phi1+j*dphi)*TMath::DegToRad();
3148 points[indx+6*n] = points[indx] = fRmin * TMath::Cos(phi);
3149 indx++;
3150 points[indx+6*n] = points[indx] = fRmin * TMath::Sin(phi);
3151 indx++;
3152 points[indx+6*n] = GetZcoord(points[indx-2], points[indx-1], dz);
3153 points[indx] = GetZcoord(points[indx-2], points[indx-1], -dz);
3154 indx++;
3155 }
3156 for (j = 0; j < n; j++) {
3157 phi = (phi1+j*dphi)*TMath::DegToRad();
3158 points[indx+6*n] = points[indx] = fRmax * TMath::Cos(phi);
3159 indx++;
3160 points[indx+6*n] = points[indx] = fRmax * TMath::Sin(phi);
3161 indx++;
3162 points[indx+6*n]= GetZcoord(points[indx-2], points[indx-1], dz);
3163 points[indx] = GetZcoord(points[indx-2], points[indx-1], -dz);
3164 indx++;
3165 }
3166 }
3167}
3168
3169////////////////////////////////////////////////////////////////////////////////
3170/// Returns numbers of vertices, segments and polygons composing the shape mesh.
3171
3172void TGeoCtub::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
3173{
3174 TGeoTubeSeg::GetMeshNumbers(nvert,nsegs,npols);
3175}
3176
3177////////////////////////////////////////////////////////////////////////////////
3178/// Return number of vertices of the mesh representation
3179
3181{
3183 Int_t numPoints = n*4;
3184 return numPoints;
3185}
3186
3187////////////////////////////////////////////////////////////////////////////////
3188/// Fills a static 3D buffer and returns a reference.
3189
3190const TBuffer3D & TGeoCtub::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
3191{
3192 static TBuffer3DCutTube buffer;
3193
3194 TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
3195
3196 if (reqSections & TBuffer3D::kShapeSpecific) {
3197 // These from TBuffer3DCutTube / TGeoCtub
3198 buffer.fRadiusInner = fRmin;
3199 buffer.fRadiusOuter = fRmax;
3200 buffer.fHalfLength = fDz;
3201 buffer.fPhiMin = fPhi1;
3202 buffer.fPhiMax = fPhi2;
3203
3204 for (UInt_t i = 0; i < 3; i++ ) {
3205 buffer.fLowPlaneNorm[i] = fNlow[i];
3206 buffer.fHighPlaneNorm[i] = fNhigh[i];
3207 }
3209 }
3210 if (reqSections & TBuffer3D::kRawSizes) {
3212 Int_t nbPnts = 4*n;
3213 Int_t nbSegs = 2*nbPnts;
3214 Int_t nbPols = nbPnts-2;
3215 if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
3217 }
3218 }
3219 if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
3220 SetPoints(buffer.fPnts);
3221 if (!buffer.fLocalFrame) {
3222 TransformPoints(buffer.fPnts, buffer.NbPnts());
3223 }
3224 SetSegsAndPols(buffer);
3226 }
3227
3228 return buffer;
3229}
3230
3231////////////////////////////////////////////////////////////////////////////////
3232/// Check the inside status for each of the points in the array.
3233/// Input: Array of point coordinates + vector size
3234/// Output: Array of Booleans for the inside of each point
3235
3236void TGeoCtub::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
3237{
3238 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
3239}
3240
3241////////////////////////////////////////////////////////////////////////////////
3242/// Compute the normal for an array o points so that norm.dot.dir is positive
3243/// Input: Arrays of point coordinates and directions + vector size
3244/// Output: Array of normal directions
3245
3246void TGeoCtub::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
3247{
3248 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
3249}
3250
3251////////////////////////////////////////////////////////////////////////////////
3252/// Compute distance from array of input points having directions specified by dirs. Store output in dists
3253
3254void TGeoCtub::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
3255{
3256 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
3257}
3258
3259////////////////////////////////////////////////////////////////////////////////
3260/// Compute distance from array of input points having directions specified by dirs. Store output in dists
3261
3262void TGeoCtub::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
3263{
3264 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
3265}
3266
3267////////////////////////////////////////////////////////////////////////////////
3268/// Compute safe distance from each of the points in the input array.
3269/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
3270/// Output: Safety values
3271
3272void TGeoCtub::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
3273{
3274 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
3275}
ROOT::R::TRInterface & r
Definition Object.C:4
#define d(i)
Definition RSha256.hxx:102
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define s1(x)
Definition RSha256.hxx:91
int Int_t
Definition RtypesCore.h:45
const Bool_t kFALSE
Definition RtypesCore.h:92
double Double_t
Definition RtypesCore.h:59
float Float_t
Definition RtypesCore.h:57
const Bool_t kTRUE
Definition RtypesCore.h:91
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:364
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition TError.cxx:187
void Warning(const char *location, const char *msgfmt,...)
Use this function in warning situations.
Definition TError.cxx:231
void Fatal(const char *location, const char *msgfmt,...)
Use this function in case of a fatal error. It will abort the program.
Definition TError.cxx:245
XFontStruct * id
Definition TGX11.cxx:109
char name[80]
Definition TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
float xmin
float ymin
float xmax
float ymax
point * points
Definition X3DBuffer.c:22
Cut tube segment description class - see TBuffer3DTypes for producer classes.
Definition TBuffer3D.h:211
Double_t fLowPlaneNorm[3]
Definition TBuffer3D.h:223
Double_t fHighPlaneNorm[3]
Definition TBuffer3D.h:224
Tube segment description class - see TBuffer3DTypes for producer classes.
Definition TBuffer3D.h:184
Double_t fPhiMax
Definition TBuffer3D.h:203
Double_t fPhiMin
Definition TBuffer3D.h:202
Complete tube description class - see TBuffer3DTypes for producer classes.
Definition TBuffer3D.h:156
Double_t fRadiusInner
Definition TBuffer3D.h:174
Double_t fRadiusOuter
Definition TBuffer3D.h:175
Double_t fHalfLength
Definition TBuffer3D.h:176
Generic 3D primitive description class.
Definition TBuffer3D.h:18
Int_t * fPols
Definition TBuffer3D.h:114
UInt_t NbPols() const
Definition TBuffer3D.h:82
UInt_t NbPnts() const
Definition TBuffer3D.h:80
UInt_t NbSegs() const
Definition TBuffer3D.h:81
Bool_t SectionsValid(UInt_t mask) const
Definition TBuffer3D.h:67
@ kShapeSpecific
Definition TBuffer3D.h:52
void SetSectionsValid(UInt_t mask)
Definition TBuffer3D.h:65
Int_t * fSegs
Definition TBuffer3D.h:113
Bool_t fLocalFrame
Definition TBuffer3D.h:90
Int_t fColor
Definition TBuffer3D.h:88
Bool_t SetRawSizes(UInt_t reqPnts, UInt_t reqPntsCapacity, UInt_t reqSegs, UInt_t reqSegsCapacity, UInt_t reqPols, UInt_t reqPolsCapacity)
Set kRaw tessellation section of buffer with supplied sizes.
Double_t * fPnts
Definition TBuffer3D.h:112
Box class.
Definition TGeoBBox.h:18
Double_t fDX
Definition TGeoBBox.h:21
virtual Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
Compute distance from outside point to surface of the box.
Definition TGeoBBox.cxx:429
virtual void InspectShape() const
Prints shape parameters.
Definition TGeoBBox.cxx:792
Double_t fOrigin[3]
Definition TGeoBBox.h:24
Double_t fDY
Definition TGeoBBox.h:22
Double_t fDZ
Definition TGeoBBox.h:23
virtual void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const
Fills the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections...
A tube segment cut with 2 planes.
Definition TGeoTube.h:169
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
computes the closest distance from given point to this shape, according to option.
virtual Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
compute distance from inside point to surface of the cut tube
virtual Bool_t Contains(const Double_t *point) const
check if point is contained in the cut tube check the lower cut plane
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
virtual void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
virtual void ComputeBBox()
compute minimum bounding box of the ctub
virtual ~TGeoCtub()
destructor
virtual void SetPoints(Double_t *points) const
Create mesh points for the cut tube.
virtual void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
in case shape has some negative parameters, these has to be computed in order to fit the mother
virtual void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
Compute safe distance from each of the points in the input array.
virtual void SetDimensions(Double_t *param)
Set dimensions of the cut tube starting from a list.
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
virtual void InspectShape() const
print shape parameters
Double_t GetZcoord(Double_t xc, Double_t yc, Double_t zc) const
compute real Z coordinate of a point belonging to either lower or higher caps (z should be either +fD...
virtual void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
Check the inside status for each of the points in the array.
virtual void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
Returns numbers of vertices, segments and polygons composing the shape mesh.
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
virtual void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
virtual Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
compute distance from outside point to surface of the cut tube
virtual Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const
Fills array with n random points located on the line segments of the shape mesh.
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
void SetCtubDimensions(Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2, Double_t lx, Double_t ly, Double_t lz, Double_t tx, Double_t ty, Double_t tz)
set dimensions of a cut tube
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide the tube along one axis.
Double_t fNlow[3]
Definition TGeoTube.h:172
Double_t fNhigh[3]
Definition TGeoTube.h:173
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Int_t GetNsegments() const
Get number of segments approximating circles.
Geometrical transformation package.
Definition TGeoMatrix.h:41
Node containing an offset.
Definition TGeoNode.h:184
Base finder class for patterns.
void SetDivIndex(Int_t index)
Base abstract class for all shapes.
Definition TGeoShape.h:26
static Double_t Big()
Definition TGeoShape.h:88
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const =0
Int_t GetBasicColor() const
Get the basic color (0-7).
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
static Double_t DistToPhiMin(const Double_t *point, const Double_t *dir, Double_t s1, Double_t c1, Double_t s2, Double_t c2, Double_t sm, Double_t cm, Bool_t in=kTRUE)
compute distance from point (inside phi) to both phi planes. Return minimum.
static Double_t SafetyPhi(const Double_t *point, Bool_t in, Double_t phi1, Double_t phi2)
Static method to compute safety w.r.t a phi corner defined by cosines/sines of the angles phi1,...
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
const char * GetPointerName() const
Provide a pointer name containing uid.
static Bool_t IsInPhiRange(const Double_t *point, Double_t phi1, Double_t phi2)
Static method to check if a point is in the phi range (phi1, phi2) [degrees].
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
virtual const char * GetName() const
Get the shape name.
static void NormalPhi(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Static method to compute normal to phi planes.
@ kGeoSavePrimitive
Definition TGeoShape.h:65
@ kGeoTubeSeg
Definition TGeoShape.h:48
@ kGeoRunTimeShape
Definition TGeoShape.h:41
static Double_t Tolerance()
Definition TGeoShape.h:91
static Bool_t IsCloseToPhi(Double_t epsil, const Double_t *point, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
True if point is closer than epsil to one of the phi planes defined by c1,s1 or c2,...
Bool_t TestShapeBit(UInt_t f) const
Definition TGeoShape.h:163
A phi segment of a tube.
Definition TGeoTube.h:89
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this tube segment shape belonging to volume "voldiv" into ndiv volumes called divname,...
virtual void Sizeof3D() const
fill size of this 3-D object
virtual void AfterStreamer()
Function called after streaming an object of this class.
virtual ~TGeoTubeSeg()
destructor
TGeoTubeSeg()
Default constructor.
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this tube segment first check if point is inside the tube
virtual void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
Double_t fPhi1
Definition TGeoTube.h:92
virtual void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
Returns numbers of vertices, segments and polygons composing the shape mesh.
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
virtual void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
virtual void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
virtual void ComputeBBox()
compute bounding box of the tube segment
Double_t fPhi2
Definition TGeoTube.h:93
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Double_t cm, Double_t sm, Double_t cdfi)
Compute distance from inside point to surface of the tube segment (static) Boundary safe algorithm.
Double_t fC2
Definition TGeoTube.h:98
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
virtual Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
compute distance from outside point to surface of the tube segment fist localize point w....
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Double_t cm, Double_t sm, Double_t cdfi)
Static method to compute distance to arbitrary tube segment from outside point Boundary safe algorith...
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Double_t fCdfi
Definition TGeoTube.h:101
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
virtual void SetPoints(Double_t *points) const
Create tube segment mesh points.
Double_t fCm
Definition TGeoTube.h:100
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2, Int_t skipz=0)
Static method to compute the closest distance from given point to this shape.
virtual Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const
Fills array with n random points located on the line segments of the shape mesh.
void InitTrigonometry()
Init frequently used trigonometric values.
void SetTubsDimensions(Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2)
Set dimensions of the tube segment.
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
Double_t fC1
Definition TGeoTube.h:96
Double_t fS1
Definition TGeoTube.h:95
Double_t fS2
Definition TGeoTube.h:97
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
computes the closest distance from given point InitTrigonometry();to this shape, according to option.
virtual void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
Check the inside status for each of the points in the array.
virtual void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
Compute safe distance from each of the points in the input array.
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons.
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
in case shape has some negative parameters, these has to be computed in order to fit the mother
virtual void InspectShape() const
print shape parameters
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Compute normal to closest surface from POINT.
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each corner
Double_t fSm
Definition TGeoTube.h:99
virtual Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
Compute distance from inside point to surface of the tube segment Boundary safe algorithm.
virtual void SetDimensions(Double_t *param)
Set dimensions of the tube segment starting from a list.
Cylindrical tube class.
Definition TGeoTube.h:18
virtual Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
Compute distance from inside point to surface of the tube Boundary safe algorithm.
Definition TGeoTube.cxx:331
virtual void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
Check the inside status for each of the points in the array.
virtual void InspectShape() const
print shape parameters
Definition TGeoTube.cxx:632
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition TGeoTube.cxx:192
virtual const char * GetAxisName(Int_t iaxis) const
Returns name of axis IAXIS.
Definition TGeoTube.cxx:543
static void DistToTube(Double_t rsq, Double_t nsq, Double_t rdotn, Double_t radius, Double_t &b, Double_t &delta)
Static method computing the distance to a tube with given radius, starting from POINT along DIR direc...
Definition TGeoTube.cxx:459
void SetTubeDimensions(Double_t rmin, Double_t rmax, Double_t dz)
Set tube dimensions.
Definition TGeoTube.cxx:909
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition TGeoTube.cxx:589
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition TGeoTube.cxx:895
virtual void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons.
Definition TGeoTube.cxx:671
Double_t fRmin
Definition TGeoTube.h:21
virtual void SetDimensions(Double_t *param)
Set tube dimensions starting from a list.
Definition TGeoTube.cxx:921
virtual Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
Compute distance from outside point to surface of the tube and safe distance Boundary safe algorithm.
Definition TGeoTube.cxx:435
Double_t fDz
Definition TGeoTube.h:23
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz)
Static method to compute distance from outside point to a tube with given parameters Boundary safe al...
Definition TGeoTube.cxx:347
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
computes the closest distance from given point to this shape, according to option.
Definition TGeoTube.cxx:826
virtual void ComputeBBox()
compute bounding box of the tube
Definition TGeoTube.cxx:209
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition TGeoTube.cxx:218
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Int_t skipz=0)
computes the closest distance from given point to this shape, according to option.
Definition TGeoTube.cxx:866
Double_t fRmax
Definition TGeoTube.h:22
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each corner
Definition TGeoTube.cxx:272
virtual Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const
Fills array with n random points located on the line segments of the shape mesh.
Definition TGeoTube.cxx:934
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
Definition TGeoTube.cxx:560
virtual ~TGeoTube()
destructor
Definition TGeoTube.cxx:185
virtual void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
virtual void Sizeof3D() const
fill size of this 3-D object
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t rmin, Double_t rmax, Double_t dz)
Compute normal to closest surface from POINT.
Definition TGeoTube.cxx:245
virtual void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
Compute safe distance from each of the points in the input array.
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this tube
Definition TGeoTube.cxx:261
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
in case shape has some negative parameters, these has to be computed in order to fit the mother
Definition TGeoTube.cxx:603
virtual void SetPoints(Double_t *points) const
create tube mesh points
Definition TGeoTube.cxx:975
virtual void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const
Compute distance from array of input points having directions specified by dirs. Store output in dist...
Bool_t HasRmin() const
Definition TGeoTube.h:69
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this tube shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition TGeoTube.cxx:481
virtual void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
Returns numbers of vertices, segments and polygons composing the shape mesh.
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
Definition TGeoTube.cxx:646
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz)
Compute distance from inside point to surface of the tube (static) Boundary safe algorithm.
Definition TGeoTube.cxx:286
TGeoTube()
Default constructor.
Definition TGeoTube.cxx:127
Volume families.
Definition TGeoVolume.h:255
void AddVolume(TGeoVolume *vol)
Add a volume with valid shape to the list of volumes.
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:49
void AddNodeOffset(TGeoVolume *vol, Int_t copy_no, Double_t offset=0, Option_t *option="")
Add a division node to the list of nodes.
TGeoMedium * GetMedium() const
Definition TGeoVolume.h:174
void SetFinder(TGeoPatternFinder *finder)
Definition TGeoVolume.h:232
Int_t GetNdaughters() const
Definition TGeoVolume.h:351
TObjArray * GetNodes()
Definition TGeoVolume.h:168
TObject * At(Int_t idx) const
Definition TObjArray.h:166
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:187
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:130
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:696
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:893
Basic string class.
Definition TString.h:136
const char * Data() const
Definition TString.h:369
return c1
Definition legend1.C:41
const Int_t n
Definition legend1.C:16
return c2
Definition legend2.C:14
Long64_t LocMin(Long64_t n, const T *a)
Return index of array with the minimum element.
Definition TMath.h:972
Short_t Max(Short_t a, Short_t b)
Definition TMathBase.h:212
T1 Sign(T1 a, T2 b)
Definition TMathBase.h:165
Double_t ATan2(Double_t y, Double_t x)
Definition TMath.h:679
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition TMath.h:1000
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition TMath.h:81
Double_t Sqrt(Double_t x)
Definition TMath.h:691
Short_t Min(Short_t a, Short_t b)
Definition TMathBase.h:180
Double_t Cos(Double_t)
Definition TMath.h:643
constexpr Double_t Pi()
Definition TMath.h:37
Double_t Sin(Double_t)
Definition TMath.h:639
constexpr Double_t RadToDeg()
Conversion from radian to degree:
Definition TMath.h:73
Short_t Abs(Short_t d)
Definition TMathBase.h:120
constexpr Double_t TwoPi()
Definition TMath.h:44
auto * t1
Definition textangle.C:20