Logo ROOT   master
Reference Guide
TGeoParaboloid.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Mihaela Gheata 20/06/04
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 TGeoParaboloid
13 \ingroup Geometry_classes
14 
15 Paraboloid class.
16 
17 A paraboloid is the solid bounded by the following surfaces:
18  - 2 planes parallel with XY cutting the Z axis at Z=-dz and Z=+dz
19  - the surface of revolution of a parabola described by:
20  `z = a*(x*x + y*y) + b`
21 
22 The parameters a and b are automatically computed from:
23  - rlo - the radius of the circle of intersection between the
24  parabolic surface and the plane z = -dz
25  - rhi - the radius of the circle of intersection between the
26  parabolic surface and the plane z = +dz
27 
28  | -dz = a*rlo*rlo + b
29  | dz = a*rhi*rhi + b where: rlo != rhi, both >= 0
30 */
31 
32 #include <iostream>
33 #include "TGeoManager.h"
34 #include "TGeoVolume.h"
35 #include "TVirtualGeoPainter.h"
36 #include "TGeoParaboloid.h"
37 #include "TBuffer3D.h"
38 #include "TBuffer3DTypes.h"
39 #include "TMath.h"
40 
42 
43 ////////////////////////////////////////////////////////////////////////////////
44 /// Dummy constructor
45 
47 {
48  fRlo = 0;
49  fRhi = 0;
50  fDz = 0;
51  fA = 0;
52  fB = 0;
54 }
55 
56 ////////////////////////////////////////////////////////////////////////////////
57 /// Default constructor specifying X and Y semiaxis length
58 
60  :TGeoBBox(0,0,0)
61 {
62  fRlo = 0;
63  fRhi = 0;
64  fDz = 0;
65  fA = 0;
66  fB = 0;
68  SetParaboloidDimensions(rlo, rhi, dz);
69  ComputeBBox();
70 }
71 
72 ////////////////////////////////////////////////////////////////////////////////
73 /// Default constructor specifying X and Y semiaxis length
74 
76  :TGeoBBox(name, 0, 0, 0)
77 {
78  fRlo = 0;
79  fRhi = 0;
80  fDz = 0;
81  fA = 0;
82  fB = 0;
84  SetParaboloidDimensions(rlo, rhi, dz);
85  ComputeBBox();
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 /// Default constructor specifying minimum and maximum radius
90 /// - param[0] = rlo
91 /// - param[1] = rhi
92 /// - param[2] = dz
93 
95 {
97  SetDimensions(param);
98  ComputeBBox();
99 }
100 
101 ////////////////////////////////////////////////////////////////////////////////
102 /// destructor
103 
105 {
106 }
107 
108 ////////////////////////////////////////////////////////////////////////////////
109 /// Computes capacity of the shape in [length^3]
110 
112 {
113  Double_t capacity = TMath::Pi()*fDz*(fRlo*fRlo+fRhi*fRhi);
114  return capacity;
115 }
116 
117 ////////////////////////////////////////////////////////////////////////////////
118 /// compute bounding box of the tube
119 
121 {
122  fDX = TMath::Max(fRlo, fRhi);
123  fDY = fDX;
124  fDZ = fDz;
125 }
126 
127 ////////////////////////////////////////////////////////////////////////////////
128 /// Compute normal to closest surface from POINT.
129 
130 void TGeoParaboloid::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
131 {
132  norm[0] = norm[1] = 0.0;
133  if (TMath::Abs(point[2]) > fDz) {
134  norm[2] = TMath::Sign(1., dir[2]);
135  return;
136  }
137  Double_t safz = fDz-TMath::Abs(point[2]);
138  Double_t r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
139  Double_t safr = TMath::Abs(r-TMath::Sqrt((point[2]-fB)/fA));
140  if (safz<safr) {
141  norm[2] = TMath::Sign(1., dir[2]);
142  return;
143  }
144  Double_t talf = -2.*fA*r;
145  Double_t calf = 1./TMath::Sqrt(1.+talf*talf);
146  Double_t salf = talf * calf;
147  Double_t phi = TMath::ATan2(point[1], point[0]);
148 
149  norm[0] = salf*TMath::Cos(phi);
150  norm[1] = salf*TMath::Sin(phi);
151  norm[2] = calf;
152  Double_t ndotd = norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2];
153  if (ndotd < 0) {
154  norm[0] = -norm[0];
155  norm[1] = -norm[1];
156  norm[2] = -norm[2];
157  }
158 }
159 
160 ////////////////////////////////////////////////////////////////////////////////
161 /// test if point is inside the elliptical tube
162 
164 {
165  if (TMath::Abs(point[2])>fDz) return kFALSE;
166  Double_t aa = fA*(point[2]-fB);
167  if (aa < 0) return kFALSE;
168  Double_t rsq = point[0]*point[0]+point[1]*point[1];
169  if (aa < fA*fA*rsq) return kFALSE;
170  return kTRUE;
171 }
172 
173 ////////////////////////////////////////////////////////////////////////////////
174 /// compute closest distance from point px,py to each vertex
175 
177 {
179  const Int_t numPoints=n*(n+1)+2;
180  return ShapeDistancetoPrimitive(numPoints, px, py);
181 }
182 
183 ////////////////////////////////////////////////////////////////////////////////
184 /// Compute distance from a point to the parabola given by:
185 /// `z = a*rsq + b; rsq = x*x+y*y`
186 
188 {
189  Double_t rsq = point[0]*point[0]+point[1]*point[1];
190  Double_t a = fA * (dir[0]*dir[0] + dir[1]*dir[1]);
191  Double_t b = 2.*fA*(point[0]*dir[0]+point[1]*dir[1])-dir[2];
192  Double_t c = fA*rsq + fB - point[2];
195  if (TMath::Abs(b)<TGeoShape::Tolerance()) return dist; // big
196  dist = -c/b;
197  if (dist < 0) return TGeoShape::Big();
198  return dist; // OK
199  }
200  Double_t ainv = 1./a;
201  Double_t sum = - b*ainv;
202  Double_t prod = c*ainv;
203  Double_t delta = sum*sum - 4.*prod;
204  if (delta<0) return dist; // big
205  delta = TMath::Sqrt(delta);
206  Double_t sone = TMath::Sign(1.,ainv);
207  Int_t i = -1;
208  while (i<2) {
209  dist = 0.5*(sum+i*sone*delta);
210  i += 2;
211  if (dist<0) continue;
212  if (dist<1.E-8) {
213  Double_t talf = -2.*fA*TMath::Sqrt(rsq);
214  Double_t phi = TMath::ATan2(point[1], point[0]);
215  Double_t ndotd = talf*(TMath::Cos(phi)*dir[0]+TMath::Sin(phi)*dir[1])+dir[2];
216  if (!in) ndotd *= -1;
217  if (ndotd<0) return dist;
218  } else return dist;
219  }
220  return TGeoShape::Big();
221 }
222 
223 ////////////////////////////////////////////////////////////////////////////////
224 /// compute distance from inside point to surface of the paraboloid
225 
226 Double_t TGeoParaboloid::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
227 {
228  if (iact<3 && safe) {
229  // compute safe distance
230  *safe = Safety(point, kTRUE);
231  if (iact==0) return TGeoShape::Big();
232  if (iact==1 && step<*safe) return TGeoShape::Big();
233  }
234 
235  Double_t dz = TGeoShape::Big();
236  if (dir[2]<0) {
237  dz = -(point[2]+fDz)/dir[2];
238  } else if (dir[2]>0) {
239  dz = (fDz-point[2])/dir[2];
240  }
241  Double_t dpara = DistToParaboloid(point, dir, kTRUE);
242  return TMath::Min(dz, dpara);
243 }
244 
245 ////////////////////////////////////////////////////////////////////////////////
246 /// compute distance from outside point to surface of the paraboloid and safe distance
247 
248 Double_t TGeoParaboloid::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
249 {
250  Double_t snxt = TGeoShape::Big();
251  if (iact<3 && safe) {
252  // compute safe distance
253  *safe = Safety(point, kFALSE);
254  if (iact==0) return TGeoShape::Big();
255  if (iact==1 && step<*safe) return TGeoShape::Big();
256  }
257  Double_t xnew, ynew, znew;
258  if (point[2]<=-fDz) {
259  if (dir[2]<=0) return TGeoShape::Big();
260  snxt = -(fDz+point[2])/dir[2];
261  // find extrapolated X and Y
262  xnew = point[0]+snxt*dir[0];
263  ynew = point[1]+snxt*dir[1];
264  if ((xnew*xnew+ynew*ynew) <= fRlo*fRlo) return snxt;
265  } else if (point[2]>=fDz) {
266  if (dir[2]>=0) return TGeoShape::Big();
267  snxt = (fDz-point[2])/dir[2];
268  // find extrapolated X and Y
269  xnew = point[0]+snxt*dir[0];
270  ynew = point[1]+snxt*dir[1];
271  if ((xnew*xnew+ynew*ynew) <= fRhi*fRhi) return snxt;
272  }
273  snxt = DistToParaboloid(point, dir, kFALSE);
274  if (snxt > 1E20) return snxt;
275  znew = point[2]+snxt*dir[2];
276  if (TMath::Abs(znew) <= fDz) return snxt;
277  return TGeoShape::Big();
278 }
279 
280 ////////////////////////////////////////////////////////////////////////////////
281 /// Divide the paraboloid along one axis.
282 
283 TGeoVolume *TGeoParaboloid::Divide(TGeoVolume * /*voldiv*/, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/,
284  Double_t /*start*/, Double_t /*step*/)
285 {
286  Error("Divide", "Paraboloid divisions not implemented");
287  return 0;
288 }
289 
290 ////////////////////////////////////////////////////////////////////////////////
291 /// Fill vector param[4] with the bounding cylinder parameters. The order
292 /// is the following : Rmin, Rmax, Phi1, Phi2
293 
295 {
296  param[0] = 0.; // Rmin
297  param[1] = fDX; // Rmax
298  param[1] *= param[1];
299  param[2] = 0.; // Phi1
300  param[3] = 360.; // Phi2
301 }
302 
303 ////////////////////////////////////////////////////////////////////////////////
304 /// in case shape has some negative parameters, these has to be computed
305 /// in order to fit the mother
306 
308 {
309  return 0;
310 }
311 
312 ////////////////////////////////////////////////////////////////////////////////
313 /// print shape parameters
314 
316 {
317  printf("*** Shape %s: TGeoParaboloid ***\n", GetName());
318  printf(" rlo = %11.5f\n", fRlo);
319  printf(" rhi = %11.5f\n", fRhi);
320  printf(" dz = %11.5f\n", fDz);
321  printf(" Bounding box:\n");
323 }
324 
325 ////////////////////////////////////////////////////////////////////////////////
326 /// Creates a TBuffer3D describing *this* shape.
327 /// Coordinates are in local reference frame.
328 
330 {
332  Int_t nbPnts = n*(n+1)+2;
333  Int_t nbSegs = n*(2*n+3);
334  Int_t nbPols = n*(n+2);
335 
337  nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 2*n*5 + n*n*6);
338 
339  if (buff)
340  {
341  SetPoints(buff->fPnts);
342  SetSegsAndPols(*buff);
343  }
344 
345  return buff;
346 }
347 
348 ////////////////////////////////////////////////////////////////////////////////
349 /// Fill TBuffer3D structure for segments and polygons.
350 
352 {
353  Int_t indx, i, j;
355 
356  Int_t c = GetBasicColor();
357 
358  Int_t nn1 = (n+1)*n+1;
359  indx = 0;
360  // Lower end-cap (n radial segments)
361  for (j=0; j<n; j++) {
362  buff.fSegs[indx++] = c+2;
363  buff.fSegs[indx++] = 0;
364  buff.fSegs[indx++] = j+1;
365  }
366  // Sectors (n)
367  for (i=0; i<n+1; i++) {
368  // lateral (circles) segments (n)
369  for (j=0; j<n; j++) {
370  buff.fSegs[indx++] = c;
371  buff.fSegs[indx++] = n*i+1+j;
372  buff.fSegs[indx++] = n*i+1+((j+1)%n);
373  }
374  if (i==n) break; // skip i=n for generators
375  // generator segments (n)
376  for (j=0; j<n; j++) {
377  buff.fSegs[indx++] = c;
378  buff.fSegs[indx++] = n*i+1+j;
379  buff.fSegs[indx++] = n*(i+1)+1+j;
380  }
381  }
382  // Upper end-cap
383  for (j=0; j<n; j++) {
384  buff.fSegs[indx++] = c+1;
385  buff.fSegs[indx++] = n*n+1+j;
386  buff.fSegs[indx++] = nn1;
387  }
388 
389  indx = 0;
390 
391  // lower end-cap (n polygons)
392  for (j=0; j<n; j++) {
393  buff.fPols[indx++] = c+2;
394  buff.fPols[indx++] = 3;
395  buff.fPols[indx++] = n+j;
396  buff.fPols[indx++] = (j+1)%n;
397  buff.fPols[indx++] = j;
398  }
399  // Sectors (n)
400  for (i=0; i<n; i++) {
401  // lateral faces (n)
402  for (j=0; j<n; j++) {
403  buff.fPols[indx++] = c;
404  buff.fPols[indx++] = 4;
405  buff.fPols[indx++] = (2*i+1)*n+j;
406  buff.fPols[indx++] = 2*(i+1)*n+j;
407  buff.fPols[indx++] = (2*i+3)*n+j;
408  buff.fPols[indx++] = 2*(i+1)*n+((j+1)%n);
409  }
410  }
411  // upper end-cap (n polygons)
412  for (j=0; j<n; j++) {
413  buff.fPols[indx++] = c+1;
414  buff.fPols[indx++] = 3;
415  buff.fPols[indx++] = 2*n*(n+1)+j;
416  buff.fPols[indx++] = 2*n*(n+1)+((j+1)%n);
417  buff.fPols[indx++] = (2*n+1)*n+j;
418  }
419 }
420 
421 ////////////////////////////////////////////////////////////////////////////////
422 /// Computes the closest distance from given point to this shape.
423 
425 {
426  Double_t safz = fDz-TMath::Abs(point[2]);
427  if (!in) safz = -safz;
428  Double_t safr = TGeoShape::Big();
429  Double_t rsq = point[0]*point[0]+point[1]*point[1];
430  Double_t z0 = fA*rsq+fB;
431  Double_t r0sq = (point[2]-fB)/fA;
432  if (r0sq<0) {
433  if (in) return 0.;
434  return safz;
435  }
436  Double_t dr = TMath::Sqrt(rsq)-TMath::Sqrt(r0sq);
437  if (in) {
438  if (dr>-1.E-8) return 0.;
439  Double_t dz = TMath::Abs(point[2]-z0);
440  safr = -dr*dz/TMath::Sqrt(dr*dr+dz*dz);
441  } else {
442  if (dr<1.E-8) return safz;
443  Double_t talf = -2.*fA*TMath::Sqrt(r0sq);
444  Double_t salf = talf/TMath::Sqrt(1.+talf*talf);
445  safr = TMath::Abs(dr*salf);
446  }
447  if (in) return TMath::Min(safr,safz);
448  return TMath::Max(safr,safz);
449 }
450 
451 ////////////////////////////////////////////////////////////////////////////////
452 /// Set paraboloid dimensions.
453 
455 {
456  if ((rlo<0) || (rhi<0) || (dz<=0) || TMath::Abs(rlo-rhi)<TGeoShape::Tolerance()) {
458  Error("SetParaboloidDimensions", "Dimensions of %s invalid: check (rlo>=0) (rhi>=0) (rlo!=rhi) dz>0",GetName());
459  return;
460  }
461  fRlo = rlo;
462  fRhi = rhi;
463  fDz = dz;
464  Double_t dd = 1./(fRhi*fRhi - fRlo*fRlo);
465  fA = 2.*fDz*dd;
466  fB = - fDz * (fRlo*fRlo + fRhi*fRhi)*dd;
467 }
468 
469 ////////////////////////////////////////////////////////////////////////////////
470 /// Set paraboloid dimensions starting from an array.
471 
473 {
474  Double_t rlo = param[0];
475  Double_t rhi = param[1];
476  Double_t dz = param[2];
477  SetParaboloidDimensions(rlo, rhi, dz);
478 }
479 
480 ////////////////////////////////////////////////////////////////////////////////
481 /// Create paraboloid mesh points.
482 /// ~~~ {.cpp}
483 /// Npoints = n*(n+1) + 2
484 /// ifirst = 0
485 /// ipoint(i,j) = 1+i*n+j; i=[0,n] j=[0,n-1]
486 /// ilast = 1+n*(n+1)
487 /// Nsegments = n*(2*n+3)
488 /// lower: (0, j+1); j=[0,n-1]
489 /// circle(i): (n*i+1+j, n*i+1+(j+1)%n); i=[0,n] j=[0,n-1]
490 /// generator(i): (n*i+1+j, n*(i+1)+1+j); i,j=[0,n-1]
491 /// upper: (n*n+1+j, (n+1)*n+1) j=[0,n-1]
492 /// Npolygons = n*(n+2)
493 /// lower: (n+j, (j+1)%n, j) j=[0,n-1]
494 /// lateral(i): ((2*i+1)*n+j, 2*(i+1)*n+j, (2*i+3)*n+j, 2*(i+1)*n+(j+1)%n)
495 /// i,j = [0,n-1]
496 /// upper: ((2n+1)*n+j, 2*n*(n+1)+(j+1)%n, 2*n*(n+1)+j) j=[0,n-1]
497 /// ~~~
498 
500 {
501  if (!points) return;
502  Double_t ttmin, ttmax;
503  ttmin = TMath::ATan2(-fDz, fRlo);
504  ttmax = TMath::ATan2(fDz, fRhi);
506  Double_t dtt = (ttmax-ttmin)/n;
507  Double_t dphi = 360./n;
508  Double_t tt;
509  Double_t r, z, delta;
510  Double_t phi, sph, cph;
511  Int_t indx = 0;
512  // center of the lower endcap:
513  points[indx++] = 0; // x
514  points[indx++] = 0; // y
515  points[indx++] = -fDz;
516  for (Int_t i=0; i<n+1; i++) { // nz planes = n+1
517  if (i==0) {
518  r = fRlo;
519  z = -fDz;
520  } else if (i==n) {
521  r = fRhi;
522  z = fDz;
523  } else {
524  tt = TMath::Tan(ttmin + i*dtt);
525  delta = tt*tt - 4*fA*fB; // should be always positive (a*b<0)
526  r = 0.5*(tt+TMath::Sqrt(delta))/fA;
527  z = r*tt;
528  }
529  for (Int_t j=0; j<n; j++) {
530  phi = j*dphi*TMath::DegToRad();
531  sph=TMath::Sin(phi);
532  cph=TMath::Cos(phi);
533  points[indx++] = r*cph;
534  points[indx++] = r*sph;
535  points[indx++] = z;
536  }
537  }
538  // center of the upper endcap
539  points[indx++] = 0; // x
540  points[indx++] = 0; // y
541  points[indx++] = fDz;
542 }
543 
544 ////////////////////////////////////////////////////////////////////////////////
545 /// Returns numbers of vertices, segments and polygons composing the shape mesh.
546 
547 void TGeoParaboloid::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
548 {
550  nvert = n*(n+1)+2;
551  nsegs = n*(2*n+3);
552  npols = n*(n+2);
553 }
554 
555 ////////////////////////////////////////////////////////////////////////////////
556 /// Returns number of vertices on the paraboloid mesh.
557 
559 {
561  return (n*(n+1)+2);
562 }
563 
564 ////////////////////////////////////////////////////////////////////////////////
565 /// Save a primitive as a C++ statement(s) on output stream "out".
566 
567 void TGeoParaboloid::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
568 {
569  if (TObject::TestBit(kGeoSavePrimitive)) return;
570  out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
571  out << " rlo = " << fRlo << ";" << std::endl;
572  out << " rhi = " << fRhi << ";" << std::endl;
573  out << " dz = " << fDZ << ";" << std::endl;
574  out << " TGeoShape *" << GetPointerName() << " = new TGeoParaboloid(\"" << GetName() << "\", rlo,rhi,dz);" << std::endl;
576 }
577 
578 ////////////////////////////////////////////////////////////////////////////////
579 /// Create paraboloid mesh points.
580 
582 {
583  if (!points) return;
584  Double_t ttmin, ttmax;
585  ttmin = TMath::ATan2(-fDz, fRlo);
586  ttmax = TMath::ATan2(fDz, fRhi);
588  Double_t dtt = (ttmax-ttmin)/n;
589  Double_t dphi = 360./n;
590  Double_t tt;
591  Double_t r, z, delta;
592  Double_t phi, sph, cph;
593  Int_t indx = 0;
594  // center of the lower endcap:
595  points[indx++] = 0; // x
596  points[indx++] = 0; // y
597  points[indx++] = -fDz;
598  for (Int_t i=0; i<n+1; i++) { // nz planes = n+1
599  if (i==0) {
600  r = fRlo;
601  z = -fDz;
602  } else if (i==n) {
603  r = fRhi;
604  z = fDz;
605  } else {
606  tt = TMath::Tan(ttmin + i*dtt);
607  delta = tt*tt - 4*fA*fB; // should be always positive (a*b<0)
608  r = 0.5*(tt+TMath::Sqrt(delta))/fA;
609  z = r*tt;
610  }
611  for (Int_t j=0; j<n; j++) {
612  phi = j*dphi*TMath::DegToRad();
613  sph=TMath::Sin(phi);
614  cph=TMath::Cos(phi);
615  points[indx++] = r*cph;
616  points[indx++] = r*sph;
617  points[indx++] = z;
618  }
619  }
620  // center of the upper endcap
621  points[indx++] = 0; // x
622  points[indx++] = 0; // y
623  points[indx++] = fDz;
624 }
625 
626 ////////////////////////////////////////////////////////////////////////////////
628 {
629 }
630 
631 ////////////////////////////////////////////////////////////////////////////////
632 /// Fills a static 3D buffer and returns a reference.
633 
634 const TBuffer3D & TGeoParaboloid::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
635 {
636  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
637  TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
638 
639  if (reqSections & TBuffer3D::kRawSizes) {
641  Int_t nbPnts = n*(n+1)+2;
642  Int_t nbSegs = n*(2*n+3);
643  Int_t nbPols = n*(n+2);
644  if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 2*n*5 + n*n*6)) {
646  }
647  }
648  if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
649  SetPoints(buffer.fPnts);
650  if (!buffer.fLocalFrame) {
651  TransformPoints(buffer.fPnts, buffer.NbPnts());
652  }
653  SetSegsAndPols(buffer);
655  }
656 
657  return buffer;
658 }
659 
660 ////////////////////////////////////////////////////////////////////////////////
661 /// Check the inside status for each of the points in the array.
662 /// Input: Array of point coordinates + vector size
663 /// Output: Array of Booleans for the inside of each point
664 
665 void TGeoParaboloid::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
666 {
667  for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
668 }
669 
670 ////////////////////////////////////////////////////////////////////////////////
671 /// Compute the normal for an array o points so that norm.dot.dir is positive
672 /// Input: Arrays of point coordinates and directions + vector size
673 /// Output: Array of normal directions
674 
675 void TGeoParaboloid::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
676 {
677  for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
678 }
679 
680 ////////////////////////////////////////////////////////////////////////////////
681 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
682 
683 void TGeoParaboloid::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
684 {
685  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
686 }
687 
688 ////////////////////////////////////////////////////////////////////////////////
689 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
690 
691 void TGeoParaboloid::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
692 {
693  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
694 }
695 
696 ////////////////////////////////////////////////////////////////////////////////
697 /// Compute safe distance from each of the points in the input array.
698 /// Input: Array of point coordinates, array of statuses for these points, size of the arrays
699 /// Output: Safety values
700 
701 void TGeoParaboloid::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
702 {
703  for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
704 }
TGeoParaboloid()
Dummy constructor.
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
static long int sum(long int i)
Definition: Factory.cxx:2272
auto * tt
Definition: textangle.C:16
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.
Box class.
Definition: TGeoBBox.h:17
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:165
virtual void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
Check the inside status for each of the points in the array.
virtual void SetPoints(Double_t *points) const
Create paraboloid mesh points.
float Float_t
Definition: RtypesCore.h:55
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
const char Option_t
Definition: RtypesCore.h:64
Geometrical transformation package.
Definition: TGeoMatrix.h:40
virtual void Sizeof3D() const
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition: TGeoVolume.h:44
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
UInt_t NbPnts() const
Definition: TBuffer3D.h:80
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoBBox.cxx:792
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
Paraboloid class.
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:695
static Double_t Tolerance()
Definition: TGeoShape.h:91
virtual void InspectShape() const
print shape parameters
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
Computes the closest distance from given point to this shape.
virtual ~TGeoParaboloid()
destructor
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide the paraboloid along one axis.
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:129
virtual void ComputeBBox()
compute bounding box of the tube
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 paraboloid and safe distance ...
Double_t fDZ
Definition: TGeoBBox.h:23
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:669
Double_t * fPnts
Definition: TBuffer3D.h:112
constexpr Double_t DegToRad()
Conversion from degree to radian: .
Definition: TMath.h:82
Double_t DistToParaboloid(const Double_t *point, const Double_t *dir, Bool_t in) const
Compute distance from a point to the parabola given by: z = a*rsq + b; rsq = x*x+y*y ...
constexpr Double_t Pi()
Definition: TMath.h:38
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
void SetSectionsValid(UInt_t mask)
Definition: TBuffer3D.h:65
virtual void SetDimensions(Double_t *param)
Set paraboloid dimensions starting from an array.
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:699
Int_t * fPols
Definition: TBuffer3D.h:114
Bool_t fLocalFrame
Definition: TBuffer3D.h:90
point * points
Definition: X3DBuffer.c:22
virtual const char * GetName() const
Get the shape name.
Definition: TGeoShape.cxx:248
virtual void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
Returns numbers of vertices, segments and polygons composing the shape mesh.
Int_t GetNsegments() const
Get number of segments approximating circles.
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
Definition: TGeoShape.cxx:259
Bool_t SectionsValid(UInt_t mask) const
Definition: TBuffer3D.h:67
Base abstract class for all shapes.
Definition: TGeoShape.h:25
ROOT::R::TRInterface & r
Definition: Object.C:4
virtual Bool_t Contains(const Double_t *point) const
test if point is inside the elliptical tube
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...
auto * a
Definition: textangle.C:12
virtual void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
Definition: TGeoShape.cxx:552
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:889
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
Generic 3D primitive description class.
Definition: TBuffer3D.h:17
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
void SetParaboloidDimensions(Double_t rlo, Double_t rhi, Double_t dz)
Set paraboloid dimensions.
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
constexpr Double_t E()
Base of natural log: .
Definition: TMath.h:97
Double_t Cos(Double_t)
Definition: TMath.h:631
const Bool_t kFALSE
Definition: RtypesCore.h:90
#define ClassImp(name)
Definition: Rtypes.h:361
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
double Double_t
Definition: RtypesCore.h:57
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
in case shape has some negative parameters, these has to be computed in order to fit the mother ...
virtual void 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...
static Double_t Big()
Definition: TGeoShape.h:88
Double_t fDY
Definition: TGeoBBox.h:22
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
Definition: TGeoShape.cxx:524
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
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 paraboloid
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:1032
virtual Int_t GetNmeshVertices() const
Returns number of vertices on the paraboloid mesh.
Int_t * fSegs
Definition: TBuffer3D.h:113
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
Double_t Sin(Double_t)
Definition: TMath.h:627
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
Double_t fDX
Definition: TGeoBBox.h:21
#define c(i)
Definition: RSha256.hxx:101
Double_t Sqrt(Double_t x)
Definition: TMath.h:681
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each vertex
const Bool_t kTRUE
Definition: RtypesCore.h:89
Int_t GetBasicColor() const
Get the basic color (0-7).
Definition: TGeoShape.cxx:673
const Int_t n
Definition: legend1.C:16
Double_t Tan(Double_t)
Definition: TMath.h:635
char name[80]
Definition: TGX11.cxx:109
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons.