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