Logo ROOT  
Reference Guide
TGeoEltu.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Mihaela Gheata 05/06/02
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/** \class TGeoEltu
13\ingroup Geometry_classes
14
15Elliptical tube class. An elliptical tube has 3 parameters
16 - A - semi-axis of the ellipse along x
17 - B - semi-axis of the ellipse along y
18 - dz - half length in z
19
20Begin_Macro(source)
21{
22 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
23 new TGeoManager("eltu", "poza6");
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->MakeEltu("ELTU",med, 30,10,40);
29 top->AddNode(vol,1);
30 gGeoManager->CloseGeometry();
31 gGeoManager->SetNsegments(50);
32 top->Draw();
33 TView *view = gPad->GetView();
34 view->ShowAxis();
35}
36End_Macro
37*/
38
39
40#include <iostream>
41
42#include "TGeoManager.h"
43#include "TGeoVolume.h"
44#include "TGeoEltu.h"
45#include "TBuffer3D.h"
46#include "TBuffer3DTypes.h"
47#include "TMath.h"
48
50
51////////////////////////////////////////////////////////////////////////////////
52/// Dummy constructor
53
55{
57}
58
59////////////////////////////////////////////////////////////////////////////////
60/// Default constructor specifying X and Y semiaxis length
61
63 :TGeoTube(0, 0, 0)
64{
66 SetEltuDimensions(a, b, dz);
68}
69
70////////////////////////////////////////////////////////////////////////////////
71/// Default constructor specifying X and Y semiaxis length
72
74 :TGeoTube(name,0.,b,dz)
75{
78 SetEltuDimensions(a, b, dz);
80}
81
82////////////////////////////////////////////////////////////////////////////////
83/// Default constructor specifying minimum and maximum radius
84/// param[0] = A
85/// param[1] = B
86/// param[2] = dz
87
89{
91 SetDimensions(param);
93}
94
95////////////////////////////////////////////////////////////////////////////////
96/// destructor
97
99{
100}
101
102////////////////////////////////////////////////////////////////////////////////
103/// Computes capacity of the shape in [length^3]
104
106{
107 Double_t capacity = 2.*TMath::Pi()*fDz*fRmin*fRmax;
108 return capacity;
109}
110
111////////////////////////////////////////////////////////////////////////////////
112/// compute bounding box of the tube
113
115{
116 fDX = fRmin;
117 fDY = fRmax;
118 fDZ = fDz;
119}
120
121////////////////////////////////////////////////////////////////////////////////
122/// Compute normal to closest surface from POINT.
123
124void TGeoEltu::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
125{
126 Double_t a = fRmin;
127 Double_t b = fRmax;
128 Double_t safr = TMath::Abs(TMath::Sqrt(point[0]*point[0]/(a*a)+point[1]*point[1]/(b*b))-1.);
129 safr *= TMath::Min(a,b);
130 Double_t safz = TMath::Abs(fDz-TMath::Abs(point[2]));
131 if (safz<safr) {
132 norm[0] = norm[1] = 0;
133 norm[2] = TMath::Sign(1.,dir[2]);
134 return;
135 }
136 norm[2] = 0.;
137 norm[0] = point[0]*b*b;
138 norm[1] = point[1]*a*a;
139 TMath::Normalize(norm);
140}
141
142////////////////////////////////////////////////////////////////////////////////
143/// test if point is inside the elliptical tube
144
146{
147 if (TMath::Abs(point[2]) > fDz) return kFALSE;
148 Double_t r2 = (point[0]*point[0])/(fRmin*fRmin)+(point[1]*point[1])/(fRmax*fRmax);
149 if (r2>1.) return kFALSE;
150 return kTRUE;
151}
152
153////////////////////////////////////////////////////////////////////////////////
154/// compute closest distance from point px,py to each vertex
155
157{
159 const Int_t numPoints=4*n;
160 return ShapeDistancetoPrimitive(numPoints, px, py);
161}
162
163////////////////////////////////////////////////////////////////////////////////
164/// compute distance from inside point to surface of the tube
165
166Double_t TGeoEltu::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
167{
170 Double_t safz1=fDz-point[2];
171 Double_t safz2=fDz+point[2];
172
173 if (iact<3 && safe) {
174 Double_t x0=TMath::Abs(point[0]);
175 Double_t y0=TMath::Abs(point[1]);
176 Double_t x1=x0;
178 Double_t y2=y0;
180 Double_t d1=(x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
181 Double_t d2=(x2-x0)*(x2-x0)+(y2-y0)*(y2-y0);
182 Double_t x3,y3;
183
184 Double_t safz = TMath::Min(safz1,safz2);
185 for (Int_t i=0; i<8; i++) {
186 if (fRmax<fRmin) {
187 x3=0.5*(x1+x2);
189 } else {
190 y3=0.5*(y1+y2);
191 x3=TMath::Sqrt((fRmax-y3)*(fRmax+y3))*fRmin/fRmax;
192 }
193 if (d1<d2) {
194 x2=x3;
195 y2=y3;
196 d2=(x2-x0)*(x2-x0)+(y2-y0)*(y2-y0);
197 } else {
198 x1=x3;
199 y1=y3;
200 d1=(x1-x0)*(x1-x0)+(y1-y0)*(y1-y0);
201 }
202 }
203 Double_t safr = TMath::Sqrt(d1)-1.0E-3;
204 *safe = TMath::Min(safz, safr);
205 if (iact==0) return TGeoShape::Big();
206 if ((iact==1) && (*safe>step)) return TGeoShape::Big();
207 }
208 // compute distance to surface
209 // Do Z
210 Double_t snxt = TGeoShape::Big();
211 if (dir[2]>0) {
212 snxt=safz1/dir[2];
213 } else {
214 if (dir[2]<0) snxt=-safz2/dir[2];
215 }
216 Double_t sz = snxt;
217 Double_t xz=point[0]+dir[0]*sz;
218 Double_t yz=point[1]+dir[1]*sz;
219 if ((xz*xz/a2+yz*yz/b2)<=1) return snxt;
220 // do elliptical surface
221 Double_t tolerance = TGeoShape::Tolerance();
222 Double_t u=dir[0]*dir[0]*b2+dir[1]*dir[1]*a2;
223 Double_t v=point[0]*dir[0]*b2+point[1]*dir[1]*a2;
224 Double_t w=point[0]*point[0]*b2+point[1]*point[1]*a2-a2*b2;
225 Double_t d=v*v-u*w;
226 if (d<0 || TGeoShape::IsSameWithinTolerance(u,0)) return tolerance;
228 snxt = (-v+sd)/u;
229
230 if (snxt<0) return tolerance;
231 return snxt;
232}
233
234////////////////////////////////////////////////////////////////////////////////
235/// compute distance from outside point to surface of the tube and safe distance
236
237Double_t TGeoEltu::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
238{
239 Double_t safz=TMath::Abs(point[2])-fDz;
242 if (iact<3 && safe) {
243 Double_t x0=TMath::Abs(point[0]);
244 Double_t y0=TMath::Abs(point[1]);
245 *safe=0.;
246 if ((x0*x0/a2+y0*y0/b2)>=1) {
247 Double_t phi1=0;
248 Double_t phi2=0.5*TMath::Pi();
249 Double_t phi3;
250 Double_t x3=0.,y3=0.,d;
251 for (Int_t i=0; i<10; i++) {
252 phi3=(phi1+phi2)*0.5;
253 x3=fRmin*TMath::Cos(phi3);
254 y3=fRmax*TMath::Sin(phi3);
255 d=y3*a2*(x0-x3)-x3*b2*(y0-y3);
256 if (d<0) phi1=phi3;
257 else phi2=phi3;
258 }
259 *safe=TMath::Sqrt((x0-x3)*(x0-x3)+(y0-y3)*(y0-y3));
260 }
261 if (safz>0) {
262 *safe=TMath::Sqrt((*safe)*(*safe)+safz*safz);
263 }
264 if (iact==0) return TGeoShape::Big();
265 if ((iact==1) && (step<*safe)) return TGeoShape::Big();
266 }
267 // compute vector distance
268 Double_t zi, tau;
269 Double_t epsil = 10.*TGeoShape::Tolerance();
270 if (safz > -epsil) {
271 // point beyond the z limit (up or down)
272 // Check if direction is outgoing
273 if (point[2]*dir[2]>0) return TGeoShape::Big();
274 // Check if direction is perpendicular to Z axis
275 if (TGeoShape::IsSameWithinTolerance(dir[2],0)) return TGeoShape::Big();
276 // select +z or -z depending on the side of the point
277 zi = (point[2] > 0) ? fDz : -fDz;
278 // Distance to zi plane position
279 tau = (zi-point[2])/dir[2];
280 // Extrapolated coordinates at the z position of the end plane.
281 Double_t xz=point[0]+dir[0]*tau;
282 Double_t yz=point[1]+dir[1]*tau;
283 if ((xz*xz/a2+yz*yz/b2)<1) return tau;
284 }
285
286// Check if the bounding box is crossed within the requested distance
287 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
288 if (sdist>=step) return TGeoShape::Big();
289 Double_t u=dir[0]*dir[0]*b2+dir[1]*dir[1]*a2; // positive
291 Double_t v=point[0]*dir[0]*b2+point[1]*dir[1]*a2;
292 Double_t w=point[0]*point[0]*b2+point[1]*point[1]*a2-a2*b2;
293 Double_t d=v*v-u*w;
294 if (d<0) return TGeoShape::Big();
296 // Biggest solution - if negative, or very close to boundary
297 // no crossing (just exiting, no re-entering possible)
298 tau = (-v+dsq)/u;
299 if (tau < epsil) return TGeoShape::Big();
300 // only entering crossing must be considered (smallest)
301 tau = (-v-dsq)/u;
302 zi=point[2]+tau*dir[2];
303 // If the crossing point is not in the Z range, there is no crossing
304 if ((TMath::Abs(zi)-fDz)>0) return TGeoShape::Big();
305 // crossing is backwards (point inside the ellipse) in Z range
306 if (tau < 0) return 0.;
307 // Point is outside and crossing the elliptical tube in Z range
308 return tau;
309}
310
311////////////////////////////////////////////////////////////////////////////////
312/// Divide the shape along one axis.
313
314TGeoVolume *TGeoEltu::Divide(TGeoVolume * /*voldiv*/, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/,
315 Double_t /*start*/, Double_t /*step*/)
316{
317 Error("Divide", "Elliptical tubes divisions not implemented");
318 return 0;
319}
320
321////////////////////////////////////////////////////////////////////////////////
322/// Fill vector param[4] with the bounding cylinder parameters. The order
323/// is the following : Rmin, Rmax, Phi1, Phi2
324
326{
327 param[0] = 0.; // Rmin
328 param[1] = TMath::Max(fRmin, fRmax); // Rmax
329 param[1] *= param[1];
330 param[2] = 0.; // Phi1
331 param[3] = 360.; // Phi2
332}
333
334////////////////////////////////////////////////////////////////////////////////
335/// in case shape has some negative parameters, these has to be computed
336/// in order to fit the mother
337
339{
340 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
341 if (!mother->TestShapeBit(kGeoEltu)) {
342 Error("GetMakeRuntimeShape", "invalid mother");
343 return 0;
344 }
345 Double_t a, b, dz;
346 a = fRmin;
347 b = fRmax;
348 dz = fDz;
349 if (fDz<0) dz=((TGeoEltu*)mother)->GetDz();
350 if (fRmin<0)
351 a = ((TGeoEltu*)mother)->GetA();
352 if (fRmax<0)
353 a = ((TGeoEltu*)mother)->GetB();
354
355 return (new TGeoEltu(a, b, dz));
356}
357
358////////////////////////////////////////////////////////////////////////////////
359/// print shape parameters
360
362{
363 printf("*** Shape %s: TGeoEltu ***\n", GetName());
364 printf(" A = %11.5f\n", fRmin);
365 printf(" B = %11.5f\n", fRmax);
366 printf(" dz = %11.5f\n", fDz);
367 printf(" Bounding box:\n");
369}
370
371////////////////////////////////////////////////////////////////////////////////
372/// computes the closest distance from given point to this shape, according
373/// to option. The matching point on the shape is stored in spoint.
374
375Double_t TGeoEltu::Safety(const Double_t *point, Bool_t /*in*/) const
376{
377 Double_t x0 = TMath::Abs(point[0]);
378 Double_t y0 = TMath::Abs(point[1]);
379 Double_t x1, y1, dx, dy;
380 Double_t safr, safz;
381 Double_t onepls = 1.+TGeoShape::Tolerance();
382 Double_t onemin = 1.-TGeoShape::Tolerance();
383 Double_t sqdist = x0*x0/(fRmin*fRmin)+y0*y0/(fRmax*fRmax);
384 Bool_t in = kTRUE;
385 if (sqdist>onepls) in = kFALSE;
386 else if (sqdist<onemin) in = kTRUE;
387 else return 0.;
388
389 if (in) {
390 x1 = fRmin*TMath::Sqrt(1.-(y0*y0)/(fRmax*fRmax));
391 y1 = fRmax*TMath::Sqrt(1.-(x0*x0)/(fRmin*fRmin));
392 dx = x1-x0;
393 dy = y1-y0;
394 if (TMath::Abs(dx)<TGeoShape::Tolerance()) return 0;
395 safr = dx*dy/TMath::Sqrt(dx*dx+dy*dy);
396 safz = fDz - TMath::Abs(point[2]);
397 return TMath::Min(safr,safz);
398 }
399
401 safr = y0 - fRmax;
402 } else {
404 safr = x0 - fRmin;
405 } else {
407 x1 = f*x0;
408 y1 = f*y0;
409 dx = x0-x1;
410 dy = y0-y1;
411 Double_t ast = fRmin*y1/fRmax;
412 Double_t bct = fRmax*x1/fRmin;
413 Double_t d = TMath::Sqrt(bct*bct+ast*ast);
414 safr = (dx*bct+dy*ast)/d;
415 }
416 }
417 safz = TMath::Abs(point[2])-fDz;
418 return TMath::Max(safr, safz);
419}
420
421////////////////////////////////////////////////////////////////////////////////
422/// Save a primitive as a C++ statement(s) on output stream "out".
423
424void TGeoEltu::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
425{
427 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
428 out << " a = " << fRmin << ";" << std::endl;
429 out << " b = " << fRmax << ";" << std::endl;
430 out << " dz = " << fDz << ";" << std::endl;
431 out << " TGeoShape *" << GetPointerName() << " = new TGeoEltu(\"" << GetName() << "\",a,b,dz);" << std::endl;
433}
434
435////////////////////////////////////////////////////////////////////////////////
436/// Set dimensions of the elliptical tube.
437
439{
440 if ((a<=0) || (b<0) || (dz<0)) {
442 }
443 fRmin=a;
444 fRmax=b;
445 fDz=dz;
446}
447
448////////////////////////////////////////////////////////////////////////////////
449/// Set shape dimensions starting from an array.
450
452{
453 Double_t a = param[0];
454 Double_t b = param[1];
455 Double_t dz = param[2];
456 SetEltuDimensions(a, b, dz);
457}
458
459////////////////////////////////////////////////////////////////////////////////
460/// Create elliptical tube mesh points
461
463{
464 Double_t dz;
465 Int_t j, n;
466
468 Double_t dphi = 360./n;
469 Double_t phi = 0;
470 Double_t cph,sph;
471 dz = fDz;
472
473 Int_t indx = 0;
474 Double_t r2,r;
477
478 if (points) {
479 for (j = 0; j < n; j++) {
480 points[indx+6*n] = points[indx] = 0;
481 indx++;
482 points[indx+6*n] = points[indx] = 0;
483 indx++;
484 points[indx+6*n] = dz;
485 points[indx] =-dz;
486 indx++;
487 }
488 for (j = 0; j < n; j++) {
489 phi = j*dphi*TMath::DegToRad();
490 sph=TMath::Sin(phi);
491 cph=TMath::Cos(phi);
492 r2=(a2*b2)/(b2+(a2-b2)*sph*sph);
493 r=TMath::Sqrt(r2);
494 points[indx+6*n] = points[indx] = r*cph;
495 indx++;
496 points[indx+6*n] = points[indx] = r*sph;
497 indx++;
498 points[indx+6*n]= dz;
499 points[indx] =-dz;
500 indx++;
501 }
502 }
503}
504
505////////////////////////////////////////////////////////////////////////////////
506/// Returns numbers of vertices, segments and polygons composing the shape mesh.
507
508void TGeoEltu::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
509{
510 TGeoTube::GetMeshNumbers(nvert,nsegs,npols);
511}
512
513////////////////////////////////////////////////////////////////////////////////
514/// Returns the number of vertices on the mesh.
515
517{
519}
520
521////////////////////////////////////////////////////////////////////////////////
522/// Create elliptical tube mesh points
523
525{
526 Double_t dz;
527 Int_t j, n;
528
530 Double_t dphi = 360./n;
531 Double_t phi = 0;
532 Double_t cph,sph;
533 dz = fDz;
534
535 Int_t indx = 0;
536 Double_t r2,r;
539
540 if (points) {
541 for (j = 0; j < n; j++) {
542 points[indx+6*n] = points[indx] = 0;
543 indx++;
544 points[indx+6*n] = points[indx] = 0;
545 indx++;
546 points[indx+6*n] = dz;
547 points[indx] =-dz;
548 indx++;
549 }
550 for (j = 0; j < n; j++) {
551 phi = j*dphi*TMath::DegToRad();
552 sph=TMath::Sin(phi);
553 cph=TMath::Cos(phi);
554 r2=(a2*b2)/(b2+(a2-b2)*sph*sph);
555 r=TMath::Sqrt(r2);
556 points[indx+6*n] = points[indx] = r*cph;
557 indx++;
558 points[indx+6*n] = points[indx] = r*sph;
559 indx++;
560 points[indx+6*n]= dz;
561 points[indx] =-dz;
562 indx++;
563 }
564 }
565}
566
567////////////////////////////////////////////////////////////////////////////////
568/// Fills a static 3D buffer and returns a reference.
569
570const TBuffer3D & TGeoEltu::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
571{
572 static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
573 TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
574
575 if (reqSections & TBuffer3D::kRawSizes) {
577 Int_t nbPnts = 4*n;
578 Int_t nbSegs = 8*n;
579 Int_t nbPols = 4*n;
580 if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
582 }
583 }
584 if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
585 SetPoints(buffer.fPnts);
586 if (!buffer.fLocalFrame) {
587 TransformPoints(buffer.fPnts, buffer.NbPnts());
588 }
589 SetSegsAndPols(buffer);
591 }
592
593 return buffer;
594}
595
596////////////////////////////////////////////////////////////////////////////////
597/// Check the inside status for each of the points in the array.
598/// Input: Array of point coordinates + vector size
599/// Output: Array of Booleans for the inside of each point
600
601void TGeoEltu::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
602{
603 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
604}
605
606////////////////////////////////////////////////////////////////////////////////
607/// Compute the normal for an array o points so that norm.dot.dir is positive
608/// Input: Arrays of point coordinates and directions + vector size
609/// Output: Array of normal directions
610
611void TGeoEltu::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
612{
613 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
614}
615
616////////////////////////////////////////////////////////////////////////////////
617/// Compute distance from array of input points having directions specified by dirs. Store output in dists
618
619void TGeoEltu::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
620{
621 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
622}
623
624////////////////////////////////////////////////////////////////////////////////
625/// Compute distance from array of input points having directions specified by dirs. Store output in dists
626
627void TGeoEltu::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
628{
629 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
630}
631
632////////////////////////////////////////////////////////////////////////////////
633/// Compute safe distance from each of the points in the input array.
634/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
635/// Output: Safety values
636
637void TGeoEltu::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
638{
639 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
640}
ROOT::R::TRInterface & r
Definition: Object.C:4
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define f(i)
Definition: RSha256.hxx:104
static const double x2[5]
static const double x1[5]
static const double x3[11]
int Int_t
Definition: RtypesCore.h:45
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
char name[80]
Definition: TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
point * points
Definition: X3DBuffer.c:22
Generic 3D primitive description class.
Definition: TBuffer3D.h:18
UInt_t NbPnts() const
Definition: TBuffer3D.h:80
Bool_t SectionsValid(UInt_t mask) const
Definition: TBuffer3D.h:67
@ kRawSizes
Definition: TBuffer3D.h:53
void SetSectionsValid(UInt_t mask)
Definition: TBuffer3D.h:65
Bool_t fLocalFrame
Definition: TBuffer3D.h:90
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
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
Elliptical tube class.
Definition: TGeoEltu.h:18
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: TGeoEltu.cxx:338
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoEltu.cxx:424
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
Definition: TGeoEltu.cxx:237
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: TGeoEltu.cxx:508
virtual Bool_t Contains(const Double_t *point) const
test if point is inside the elliptical tube
Definition: TGeoEltu.cxx:145
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoEltu.cxx:124
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoEltu.cxx:325
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
Definition: TGeoEltu.cxx:570
virtual ~TGeoEltu()
destructor
Definition: TGeoEltu.cxx:98
virtual void SetDimensions(Double_t *param)
Set shape dimensions starting from an array.
Definition: TGeoEltu.cxx:451
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: TGeoEltu.cxx:601
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: TGeoEltu.cxx:627
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: TGeoEltu.cxx:619
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
Definition: TGeoEltu.cxx:166
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: TGeoEltu.cxx:375
TGeoEltu()
Dummy constructor.
Definition: TGeoEltu.cxx:54
virtual void InspectShape() const
print shape parameters
Definition: TGeoEltu.cxx:361
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each vertex
Definition: TGeoEltu.cxx:156
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: TGeoEltu.cxx:637
virtual void ComputeBBox()
compute bounding box of the tube
Definition: TGeoEltu.cxx:114
virtual Int_t GetNmeshVertices() const
Returns the number of vertices on the mesh.
Definition: TGeoEltu.cxx:516
void SetEltuDimensions(Double_t a, Double_t b, Double_t dz)
Set dimensions of the elliptical tube.
Definition: TGeoEltu.cxx:438
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoEltu.cxx:105
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide the shape along one axis.
Definition: TGeoEltu.cxx:314
virtual void SetPoints(Double_t *points) const
Create elliptical tube mesh points.
Definition: TGeoEltu.cxx:462
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: TGeoEltu.cxx:611
Int_t GetNsegments() const
Get number of segments approximating circles.
Geometrical transformation package.
Definition: TGeoMatrix.h:41
Base abstract class for all shapes.
Definition: TGeoShape.h:26
static Double_t Big()
Definition: TGeoShape.h:88
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 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
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
@ kGeoSavePrimitive
Definition: TGeoShape.h:65
@ kGeoEltu
Definition: TGeoShape.h:54
@ kGeoRunTimeShape
Definition: TGeoShape.h:41
static Double_t Tolerance()
Definition: TGeoShape.h:91
Bool_t TestShapeBit(UInt_t f) const
Definition: TGeoShape.h:163
Cylindrical tube class.
Definition: TGeoTube.h:18
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons.
Definition: TGeoTube.cxx:669
Double_t fRmin
Definition: TGeoTube.h:21
Double_t fDz
Definition: TGeoTube.h:23
Double_t fRmax
Definition: TGeoTube.h:22
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
Definition: TGeoTube.cxx:1099
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
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition: TGeoVolume.h:49
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
const Int_t n
Definition: legend1.C:16
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
Float_t Normalize(Float_t v[3])
Normalize a vector v in place.
Definition: TMath.cxx:495
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:165
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
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
auto * a
Definition: textangle.C:12