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