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