Logo ROOT  
Reference Guide
TGeoCone.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 31/01/02
3 // TGeoCone::Contains() and DistFromInside() implemented by Mihaela Gheata
4 
5 /*************************************************************************
6  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7  * All rights reserved. *
8  * *
9  * For the licensing terms see $ROOTSYS/LICENSE. *
10  * For the list of contributors see $ROOTSYS/README/CREDITS. *
11  *************************************************************************/
12 
13 /** \class TGeoCone
14 \ingroup Geometry_classes
15 
16 Conical tube class. It has 5 parameters :
17  - dz - half length in z
18  - Rmin1, Rmax1 - inside and outside radii at -dz
19  - Rmin2, Rmax2 - inside and outside radii at +dz
20 
21 Begin_Macro(source)
22 {
23  TCanvas *c = new TCanvas("c", "c",0,0,600,600);
24  new TGeoManager("cone", "poza4");
25  TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
26  TGeoMedium *med = new TGeoMedium("MED",1,mat);
27  TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
28  gGeoManager->SetTopVolume(top);
29  TGeoVolume *vol = gGeoManager->MakeCone("CONE",med, 40,10,20,35,45);
30  vol->SetLineWidth(2);
31  top->AddNode(vol,1);
32  gGeoManager->CloseGeometry();
33  gGeoManager->SetNsegments(30);
34  top->Draw();
35  TView *view = gPad->GetView();
36  view->ShowAxis();
37 }
38 End_Macro
39 */
40 
41 
42 /** \class TGeoConeSeg
43 \ingroup Geometry_classes
44 
45 A phi segment of a conical tube. Has 7 parameters :
46  - the same 5 as a cone;
47  - first phi limit (in degrees)
48  - second phi limit
49 
50 Begin_Macro(source)
51 {
52  TCanvas *c = new TCanvas("c", "c",0,0,600,600);
53  new TGeoManager("coneseg", "poza5");
54  TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
55  TGeoMedium *med = new TGeoMedium("MED",1,mat);
56  TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
57  gGeoManager->SetTopVolume(top);
58  TGeoVolume *vol = gGeoManager->MakeCons("CONESEG",med, 40,30,40,10,20,-30,250);
59  top->AddNode(vol,1);
60  gGeoManager->CloseGeometry();
61  gGeoManager->SetNsegments(30);
62  top->Draw();
63  TView *view = gPad->GetView();
64  view->ShowAxis();
65 }
66 End_Macro
67 */
68 
69 #include <iostream>
70 
71 #include "TGeoManager.h"
72 #include "TGeoVolume.h"
73 #include "TVirtualGeoPainter.h"
74 #include "TGeoCone.h"
75 #include "TBuffer3D.h"
76 #include "TBuffer3DTypes.h"
77 #include "TMath.h"
78 
80 
81 ////////////////////////////////////////////////////////////////////////////////
82 /// Default constructor
83 
85 {
87  fDz = 0.0;
88  fRmin1 = 0.0;
89  fRmax1 = 0.0;
90  fRmin2 = 0.0;
91  fRmax2 = 0.0;
92 }
93 
94 ////////////////////////////////////////////////////////////////////////////////
95 /// Default constructor specifying minimum and maximum radius
96 
98  Double_t rmin2, Double_t rmax2)
99  :TGeoBBox(0, 0, 0)
100 {
102  SetConeDimensions(dz, rmin1, rmax1, rmin2, rmax2);
103  if ((dz<0) || (rmin1<0) || (rmax1<0) || (rmin2<0) || (rmax2<0)) {
105  }
106  else ComputeBBox();
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 /// Default constructor specifying minimum and maximum radius
111 
112 TGeoCone::TGeoCone(const char *name, Double_t dz, Double_t rmin1, Double_t rmax1,
113  Double_t rmin2, Double_t rmax2)
114  :TGeoBBox(name, 0, 0, 0)
115 {
117  SetConeDimensions(dz, rmin1, rmax1, rmin2, rmax2);
118  if ((dz<0) || (rmin1<0) || (rmax1<0) || (rmin2<0) || (rmax2<0)) {
120  }
121  else ComputeBBox();
122 }
123 
124 ////////////////////////////////////////////////////////////////////////////////
125 /// Default constructor specifying minimum and maximum radius
126 /// - param[0] = dz
127 /// - param[1] = Rmin1
128 /// - param[2] = Rmax1
129 /// - param[3] = Rmin2
130 /// - param[4] = Rmax2
131 
133  :TGeoBBox(0, 0, 0)
134 {
136  SetDimensions(param);
137  if ((fDz<0) || (fRmin1<0) || (fRmax1<0) || (fRmin2<0) || (fRmax2<0))
139  else ComputeBBox();
140 }
141 
142 ////////////////////////////////////////////////////////////////////////////////
143 /// Computes capacity of the shape in [length^3]
144 
146 {
148 }
149 
150 ////////////////////////////////////////////////////////////////////////////////
151 /// Computes capacity of the shape in [length^3]
152 
154 {
155  Double_t capacity = (2.*dz*TMath::Pi()/3.)*(rmax1*rmax1+rmax2*rmax2+rmax1*rmax2-
156  rmin1*rmin1-rmin2*rmin2-rmin1*rmin2);
157  return capacity;
158 }
159 
160 ////////////////////////////////////////////////////////////////////////////////
161 /// destructor
162 
164 {
165 }
166 
167 ////////////////////////////////////////////////////////////////////////////////
168 /// compute bounding box of the sphere
169 
171 {
172  TGeoBBox *box = (TGeoBBox*)this;
173  box->SetBoxDimensions(TMath::Max(fRmax1, fRmax2), TMath::Max(fRmax1, fRmax2), fDz);
174  memset(fOrigin, 0, 3*sizeof(Double_t));
175 }
176 
177 ////////////////////////////////////////////////////////////////////////////////
178 /// Compute normal to closest surface from POINT.
179 
180 void TGeoCone::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
181 {
182  Double_t safr,safe,phi;
183  memset(norm,0,3*sizeof(Double_t));
184  phi = TMath::ATan2(point[1],point[0]);
185  Double_t cphi = TMath::Cos(phi);
186  Double_t sphi = TMath::Sin(phi);
187  Double_t ro1 = 0.5*(fRmin1+fRmin2);
188  Double_t tg1 = 0.5*(fRmin2-fRmin1)/fDz;
189  Double_t cr1 = 1./TMath::Sqrt(1.+tg1*tg1);
190  Double_t ro2 = 0.5*(fRmax1+fRmax2);
191  Double_t tg2 = 0.5*(fRmax2-fRmax1)/fDz;
192  Double_t cr2 = 1./TMath::Sqrt(1.+tg2*tg2);
193 
194  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
195  Double_t rin = tg1*point[2]+ro1;
196  Double_t rout = tg2*point[2]+ro2;
197  safe = TMath::Abs(fDz-TMath::Abs(point[2]));
198  norm[2] = 1;
199 
200  safr = (ro1>0)?(TMath::Abs((r-rin)*cr1)):TGeoShape::Big();
201  if (safr<safe) {
202  safe = safr;
203  norm[0] = cr1*cphi;
204  norm[1] = cr1*sphi;
205  norm[2] = -tg1*cr1;
206  }
207  safr = TMath::Abs((rout-r)*cr2);
208  if (safr<safe) {
209  norm[0] = cr2*cphi;
210  norm[1] = cr2*sphi;
211  norm[2] = -tg2*cr2;
212  }
213  if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
214  norm[0] = -norm[0];
215  norm[1] = -norm[1];
216  norm[2] = -norm[2];
217  }
218 }
219 
220 ////////////////////////////////////////////////////////////////////////////////
221 /// Compute normal to closest surface from POINT.
222 
223 void TGeoCone::ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm,
224  Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
225 {
226  Double_t safe,phi;
227  memset(norm,0,3*sizeof(Double_t));
228  phi = TMath::ATan2(point[1],point[0]);
229  Double_t cphi = TMath::Cos(phi);
230  Double_t sphi = TMath::Sin(phi);
231  Double_t ro1 = 0.5*(rmin1+rmin2);
232  Double_t tg1 = 0.5*(rmin2-rmin1)/dz;
233  Double_t cr1 = 1./TMath::Sqrt(1.+tg1*tg1);
234  Double_t ro2 = 0.5*(rmax1+rmax2);
235  Double_t tg2 = 0.5*(rmax2-rmax1)/dz;
236  Double_t cr2 = 1./TMath::Sqrt(1.+tg2*tg2);
237 
238  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
239  Double_t rin = tg1*point[2]+ro1;
240  Double_t rout = tg2*point[2]+ro2;
241  safe = (ro1>0)?(TMath::Abs((r-rin)*cr1)):TGeoShape::Big();
242  norm[0] = cr1*cphi;
243  norm[1] = cr1*sphi;
244  norm[2] = -tg1*cr1;
245  if (TMath::Abs((rout-r)*cr2)<safe) {
246  norm[0] = cr2*cphi;
247  norm[1] = cr2*sphi;
248  norm[2] = -tg2*cr2;
249  }
250  if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
251  norm[0] = -norm[0];
252  norm[1] = -norm[1];
253  norm[2] = -norm[2];
254  }
255 }
256 
257 ////////////////////////////////////////////////////////////////////////////////
258 /// test if point is inside this cone
259 
260 Bool_t TGeoCone::Contains(const Double_t *point) const
261 {
262  if (TMath::Abs(point[2]) > fDz) return kFALSE;
263  Double_t r2 = point[0]*point[0]+point[1]*point[1];
264  Double_t rl = 0.5*(fRmin2*(point[2]+fDz)+fRmin1*(fDz-point[2]))/fDz;
265  Double_t rh = 0.5*(fRmax2*(point[2]+fDz)+fRmax1*(fDz-point[2]))/fDz;
266  if ((r2<rl*rl) || (r2>rh*rh)) return kFALSE;
267  return kTRUE;
268 }
269 
270 ////////////////////////////////////////////////////////////////////////////////
271 /// Compute distance from inside point to surface of the cone (static)
272 /// Boundary safe algorithm.
273 
275  Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
276 {
277  if (dz<=0) return TGeoShape::Big();
278  // compute distance to surface
279  // Do Z
280  Double_t sz = TGeoShape::Big();
281  if (dir[2]) {
282  sz = (TMath::Sign(dz, dir[2])-point[2])/dir[2];
283  if (sz<=0) return 0.0;
284  }
285  Double_t rsq=point[0]*point[0]+point[1]*point[1];
286  Double_t zinv = 1./dz;
287  Double_t rin = 0.5*(rmin1+rmin2+(rmin2-rmin1)*point[2]*zinv);
288  // Do Rmin
290  Double_t b,delta,zi;
291  if (rin>0) {
292  // Protection in case point is actually outside the cone
293  if (rsq < rin*(rin+TGeoShape::Tolerance())) {
294  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]+0.5*(rmin1-rmin2)*dir[2]*zinv*TMath::Sqrt(rsq);
295  if (ddotn<=0) return 0.0;
296  } else {
297  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
298  if (delta>0) {
299  sr = -b-delta;
300  if (sr>0) {
301  zi = point[2]+sr*dir[2];
302  if (TMath::Abs(zi)<=dz) return TMath::Min(sz,sr);
303  }
304  sr = -b+delta;
305  if (sr>0) {
306  zi = point[2]+sr*dir[2];
307  if (TMath::Abs(zi)<=dz) return TMath::Min(sz,sr);
308  }
309  }
310  }
311  }
312  // Do Rmax
313  Double_t rout = 0.5*(rmax1+rmax2+(rmax2-rmax1)*point[2]*zinv);
314  if (rsq > rout*(rout-TGeoShape::Tolerance())) {
315  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]+0.5*(rmax1-rmax2)*dir[2]*zinv*TMath::Sqrt(rsq);
316  if (ddotn>=0) return 0.0;
317  TGeoCone::DistToCone(point, dir, dz, rmax1, rmax2, b, delta);
318  if (delta<0) return 0.0;
319  sr = -b+delta;
320  if (sr<0) return sz;
321  if (TMath::Abs(-b-delta)>sr) return sz;
322  zi = point[2]+sr*dir[2];
323  if (TMath::Abs(zi)<=dz) return TMath::Min(sz,sr);
324  return sz;
325  }
326  TGeoCone::DistToCone(point, dir, dz, rmax1, rmax2, b, delta);
327  if (delta>0) {
328  sr = -b-delta;
329  if (sr>0) {
330  zi = point[2]+sr*dir[2];
331  if (TMath::Abs(zi)<=dz) return TMath::Min(sz,sr);
332  }
333  sr = -b+delta;
334  if (sr>TGeoShape::Tolerance()) {
335  zi = point[2]+sr*dir[2];
336  if (TMath::Abs(zi)<=dz) return TMath::Min(sz,sr);
337  }
338  }
339  return sz;
340 }
341 
342 ////////////////////////////////////////////////////////////////////////////////
343 /// Compute distance from inside point to surface of the cone
344 /// Boundary safe algorithm.
345 
346 Double_t TGeoCone::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
347 {
348  if (iact<3 && safe) {
349  *safe = Safety(point, kTRUE);
350  if (iact==0) return TGeoShape::Big();
351  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
352  }
353  // compute distance to surface
354  return TGeoCone::DistFromInsideS(point, dir, fDz, fRmin1, fRmax1, fRmin2, fRmax2);
355 }
356 
357 ////////////////////////////////////////////////////////////////////////////////
358 /// Compute distance from outside point to surface of the tube
359 /// Boundary safe algorithm.
360 
362  Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
363 {
364  // compute distance to Z planes
365  if (dz<=0) return TGeoShape::Big();
366  Double_t snxt;
367  Double_t xp, yp, zp;
368  Bool_t inz = kTRUE;
369 
370  if (point[2]<=-dz) {
371  if (dir[2]<=0) return TGeoShape::Big();
372  snxt = (-dz-point[2])/dir[2];
373  xp = point[0]+snxt*dir[0];
374  yp = point[1]+snxt*dir[1];
375  Double_t r2 = xp*xp+yp*yp;
376  if ((r2>=rmin1*rmin1) && (r2<=rmax1*rmax1)) return snxt;
377  inz = kFALSE;
378  } else {
379  if (point[2]>=dz) {
380  if (dir[2]>=0) return TGeoShape::Big();
381  snxt = (dz-point[2])/dir[2];
382  xp = point[0]+snxt*dir[0];
383  yp = point[1]+snxt*dir[1];
384  Double_t r2 = xp*xp+yp*yp;
385  if ((r2>=rmin2*rmin2) && (r2<=rmax2*rmax2)) return snxt;
386  inz = kFALSE;
387  }
388  }
389 
390  Double_t rsq = point[0]*point[0]+point[1]*point[1];
391  Double_t dzinv = 1./dz;
392  Double_t ro1=0.5*(rmin1+rmin2);
393  Bool_t hasrmin = (ro1>0)?kTRUE:kFALSE;
394  Double_t tg1 = 0.;
395  Double_t rin = 0.;
396  Bool_t inrmin = kTRUE; // r>=rmin
397  if (hasrmin) {
398  tg1=0.5*(rmin2-rmin1)*dzinv;
399  rin=ro1+tg1*point[2];
400  if (rin>0 && rsq<rin*(rin-TGeoShape::Tolerance())) inrmin=kFALSE;
401  }
402  Double_t ro2=0.5*(rmax1+rmax2);
403  Double_t tg2=0.5*(rmax2-rmax1)*dzinv;
404  Double_t rout=tg2*point[2]+ro2;
405  Bool_t inrmax = kFALSE;
406  if (rout>0 && rsq<rout*(rout+TGeoShape::Tolerance())) inrmax=kTRUE;
407  Bool_t in = inz & inrmin & inrmax;
408  Double_t b,delta;
409  // If inside cone, we are most likely on a boundary within machine precision.
410  if (in) {
411  Double_t r=TMath::Sqrt(rsq);
412  Double_t safz = dz-TMath::Abs(point[2]); // positive
413  Double_t safrmin = (hasrmin)?(r-rin):TGeoShape::Big();
414  Double_t safrmax = rout - r;
415  if (safz<=safrmin && safz<=safrmax) {
416  // on Z boundary
417  if (point[2]*dir[2]<0) return 0.0;
418  return TGeoShape::Big();
419  }
420  if (safrmax<safrmin) {
421  // on rmax boundary
422  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]-tg2*dir[2]*r;
423  if (ddotn<=0) return 0.0;
424  return TGeoShape::Big();
425  }
426  // on rmin boundary
427  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]-tg1*dir[2]*r;
428  if (ddotn>=0) return 0.0;
429  // we can cross (+) solution of rmin
430  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
431 
432  if (delta<0) return 0.0;
433  snxt = -b+delta;
434  if (snxt<0) return TGeoShape::Big();
435  if (TMath::Abs(-b-delta)>snxt) return TGeoShape::Big();
436  zp = point[2]+snxt*dir[2];
437  if (TMath::Abs(zp)<=dz) return snxt;
438  return TGeoShape::Big();
439  }
440 
441  // compute distance to inner cone
442  snxt = TGeoShape::Big();
443  if (!inrmin) {
444  // ray can cross inner cone (but not only!)
445  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
446  if (delta<0) return TGeoShape::Big();
447  snxt = -b+delta;
448  if (snxt>0) {
449  zp = point[2]+snxt*dir[2];
450  if (TMath::Abs(zp)<=dz) return snxt;
451  }
452  snxt = -b-delta;
453  if (snxt>0) {
454  zp = point[2]+snxt*dir[2];
455  if (TMath::Abs(zp)<=dz) return snxt;
456  }
457  snxt = TGeoShape::Big();
458  } else {
459  if (hasrmin) {
460  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
461  if (delta>0) {
462  Double_t din = -b+delta;
463  if (din>0) {
464  zp = point[2]+din*dir[2];
465  if (TMath::Abs(zp)<=dz) snxt = din;
466  }
467  }
468  }
469  }
470 
471  if (inrmax) return snxt;
472  // We can cross outer cone, both solutions possible
473  // compute distance to outer cone
474  TGeoCone::DistToCone(point, dir, dz, rmax1, rmax2, b, delta);
475  if (delta<0) return snxt;
476  Double_t dout = -b-delta;
477  if (dout>0 && dout<snxt) {
478  zp = point[2]+dout*dir[2];
479  if (TMath::Abs(zp)<=dz) return dout;
480  }
481  dout = -b+delta;
482  if (dout<=0 || dout>snxt) return snxt;
483  zp = point[2]+dout*dir[2];
484  if (TMath::Abs(zp)<=dz) return dout;
485  return snxt;
486 }
487 
488 ////////////////////////////////////////////////////////////////////////////////
489 /// compute distance from outside point to surface of the tube
490 
491 Double_t TGeoCone::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
492 {
493  // compute safe radius
494  if (iact<3 && safe) {
495  *safe = Safety(point, kFALSE);
496  if (iact==0) return TGeoShape::Big();
497  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
498  }
499  // Check if the bounding box is crossed within the requested distance
500  Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
501  if (sdist>=step) return TGeoShape::Big();
502  // compute distance to Z planes
503  return TGeoCone::DistFromOutsideS(point, dir, fDz, fRmin1, fRmax1, fRmin2, fRmax2);
504 }
505 
506 ////////////////////////////////////////////////////////////////////////////////
507 /// Static method to compute distance to a conical surface with :
508 /// - r1, z1 - radius and Z position of lower base
509 /// - r2, z2 - radius and Z position of upper base
510 
511 void TGeoCone::DistToCone(const Double_t *point, const Double_t *dir, Double_t dz, Double_t r1, Double_t r2,
512  Double_t &b, Double_t &delta)
513 {
514  delta = -1.;
515  if (dz<0) return;
516  Double_t ro0 = 0.5*(r1+r2);
517  Double_t tz = 0.5*(r2-r1)/dz;
518  Double_t rsq = point[0]*point[0] + point[1]*point[1];
519  Double_t rc = ro0 + point[2]*tz;
520 
521  Double_t a = dir[0]*dir[0] + dir[1]*dir[1] - tz*tz*dir[2]*dir[2];
522  b = point[0]*dir[0] + point[1]*dir[1] - tz*rc*dir[2];
523  Double_t c = rsq - rc*rc;
524 
526  if (TMath::Abs(b)<TGeoShape::Tolerance()) return;
527  b = 0.5*c/b;
528  delta = 0.;
529  return;
530  }
531  a = 1./a;
532  b *= a;
533  c *= a;
534  delta = b*b - c;
535  if (delta>0) {
536  delta = TMath::Sqrt(delta);
537  } else {
538  delta = -1.;
539  }
540 }
541 
542 ////////////////////////////////////////////////////////////////////////////////
543 /// compute closest distance from point px,py to each corner
544 
546 {
548  const Int_t numPoints = 4*n;
549  return ShapeDistancetoPrimitive(numPoints, px, py);
550 }
551 
552 ////////////////////////////////////////////////////////////////////////////////
553 /// Divide this cone shape belonging to volume "voldiv" into ndiv volumes
554 /// called divname, from start position with the given step. Returns pointer
555 /// to created division cell volume in case of Z divisions. For Z division
556 /// creates all volumes with different shapes and returns pointer to volume that
557 /// was divided. In case a wrong division axis is supplied, returns pointer to
558 /// volume that was divided.
559 
560 TGeoVolume *TGeoCone::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
561  Double_t start, Double_t step)
562 {
563  TGeoShape *shape; //--- shape to be created
564  TGeoVolume *vol; //--- division volume to be created
565  TGeoVolumeMulti *vmulti; //--- generic divided volume
566  TGeoPatternFinder *finder; //--- finder to be attached
567  TString opt = ""; //--- option to be attached
568  Int_t id;
569  Double_t end = start+ndiv*step;
570  switch (iaxis) {
571  case 1: //--- R division
572  Error("Divide","division of a cone on R not implemented");
573  return 0;
574  case 2: // --- Phi division
575  finder = new TGeoPatternCylPhi(voldiv, ndiv, start, end);
576  voldiv->SetFinder(finder);
577  finder->SetDivIndex(voldiv->GetNdaughters());
578  shape = new TGeoConeSeg(fDz, fRmin1, fRmax1, fRmin2, fRmax2, -step/2, step/2);
579  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
580  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
581  vmulti->AddVolume(vol);
582  opt = "Phi";
583  for (id=0; id<ndiv; id++) {
584  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
585  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
586  }
587  return vmulti;
588  case 3: //--- Z division
589  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
590  finder = new TGeoPatternZ(voldiv, ndiv, start, end);
591  voldiv->SetFinder(finder);
592  finder->SetDivIndex(voldiv->GetNdaughters());
593  for (id=0; id<ndiv; id++) {
594  Double_t z1 = start+id*step;
595  Double_t z2 = start+(id+1)*step;
596  Double_t rmin1n = 0.5*(fRmin1*(fDz-z1)+fRmin2*(fDz+z1))/fDz;
597  Double_t rmax1n = 0.5*(fRmax1*(fDz-z1)+fRmax2*(fDz+z1))/fDz;
598  Double_t rmin2n = 0.5*(fRmin1*(fDz-z2)+fRmin2*(fDz+z2))/fDz;
599  Double_t rmax2n = 0.5*(fRmax1*(fDz-z2)+fRmax2*(fDz+z2))/fDz;
600  shape = new TGeoCone(0.5*step,rmin1n, rmax1n, rmin2n, rmax2n);
601  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
602  vmulti->AddVolume(vol);
603  opt = "Z";
604  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
605  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
606  }
607  return vmulti;
608  default:
609  Error("Divide", "Wrong axis type for division");
610  return 0;
611  }
612 }
613 
614 ////////////////////////////////////////////////////////////////////////////////
615 /// Returns name of axis IAXIS.
616 
617 const char *TGeoCone::GetAxisName(Int_t iaxis) const
618 {
619  switch (iaxis) {
620  case 1:
621  return "R";
622  case 2:
623  return "PHI";
624  case 3:
625  return "Z";
626  default:
627  return "undefined";
628  }
629 }
630 
631 ////////////////////////////////////////////////////////////////////////////////
632 /// Get range of shape for a given axis.
633 
635 {
636  xlo = 0;
637  xhi = 0;
638  Double_t dx = 0;
639  switch (iaxis) {
640  case 2:
641  xlo = 0.;
642  xhi = 360.;
643  return 360.;
644  case 3:
645  xlo = -fDz;
646  xhi = fDz;
647  dx = xhi-xlo;
648  return dx;
649  }
650  return dx;
651 }
652 
653 ////////////////////////////////////////////////////////////////////////////////
654 /// Fill vector param[4] with the bounding cylinder parameters. The order
655 /// is the following : Rmin, Rmax, Phi1, Phi2, dZ
656 
658 {
659  param[0] = TMath::Min(fRmin1, fRmin2); // Rmin
660  param[0] *= param[0];
661  param[1] = TMath::Max(fRmax1, fRmax2); // Rmax
662  param[1] *= param[1];
663  param[2] = 0.; // Phi1
664  param[3] = 360.; // Phi2
665 }
666 
667 ////////////////////////////////////////////////////////////////////////////////
668 /// in case shape has some negative parameters, these has to be computed
669 /// in order to fit the mother
670 
672 {
673  if (!TestShapeBit(kGeoRunTimeShape)) return 0;
674  if (!mother->TestShapeBit(kGeoCone)) {
675  Error("GetMakeRuntimeShape", "invalid mother");
676  return 0;
677  }
678  Double_t rmin1, rmax1, rmin2, rmax2, dz;
679  rmin1 = fRmin1;
680  rmax1 = fRmax1;
681  rmin2 = fRmin2;
682  rmax2 = fRmax2;
683  dz = fDz;
684  if (fDz<0) dz=((TGeoCone*)mother)->GetDz();
685  if (fRmin1<0)
686  rmin1 = ((TGeoCone*)mother)->GetRmin1();
687  if (fRmax1<0)
688  rmax1 = ((TGeoCone*)mother)->GetRmax1();
689  if (fRmin2<0)
690  rmin2 = ((TGeoCone*)mother)->GetRmin2();
691  if (fRmax2<0)
692  rmax2 = ((TGeoCone*)mother)->GetRmax2();
693 
694  return (new TGeoCone(GetName(), dz, rmin1, rmax1, rmin2, rmax2));
695 }
696 
697 ////////////////////////////////////////////////////////////////////////////////
698 /// Fills array with n random points located on the line segments of the shape mesh.
699 /// The output array must be provided with a length of minimum 3*npoints. Returns
700 /// true if operation is implemented.
701 
703 {
704  if (npoints > (npoints/2)*2) {
705  Error("GetPointsOnSegments","Npoints must be even number");
706  return kFALSE;
707  }
708  Bool_t hasrmin = (fRmin1>0 || fRmin2>0)?kTRUE:kFALSE;
709  Int_t nc = 0;
710  if (hasrmin) nc = (Int_t)TMath::Sqrt(0.5*npoints);
711  else nc = (Int_t)TMath::Sqrt(1.*npoints);
712  Double_t dphi = TMath::TwoPi()/nc;
713  Double_t phi = 0;
714  Int_t ntop = 0;
715  if (hasrmin) ntop = npoints/2 - nc*(nc-1);
716  else ntop = npoints - nc*(nc-1);
717  Double_t dz = 2*fDz/(nc-1);
718  Double_t z = 0;
719  Int_t icrt = 0;
720  Int_t nphi = nc;
721  Double_t rmin = 0.;
722  Double_t rmax = 0.;
723  // loop z sections
724  for (Int_t i=0; i<nc; i++) {
725  if (i == (nc-1)) nphi = ntop;
726  z = -fDz + i*dz;
727  if (hasrmin) rmin = 0.5*(fRmin1+fRmin2) + 0.5*(fRmin2-fRmin1)*z/fDz;
728  rmax = 0.5*(fRmax1+fRmax2) + 0.5*(fRmax2-fRmax1)*z/fDz;
729  // loop points on circle sections
730  for (Int_t j=0; j<nphi; j++) {
731  phi = j*dphi;
732  if (hasrmin) {
733  array[icrt++] = rmin * TMath::Cos(phi);
734  array[icrt++] = rmin * TMath::Sin(phi);
735  array[icrt++] = z;
736  }
737  array[icrt++] = rmax * TMath::Cos(phi);
738  array[icrt++] = rmax * TMath::Sin(phi);
739  array[icrt++] = z;
740  }
741  }
742  return kTRUE;
743 }
744 
745 
746 ////////////////////////////////////////////////////////////////////////////////
747 /// print shape parameters
748 
750 {
751  printf("*** Shape %s TGeoCone ***\n", GetName());
752  printf(" dz =: %11.5f\n", fDz);
753  printf(" Rmin1 = %11.5f\n", fRmin1);
754  printf(" Rmax1 = %11.5f\n", fRmax1);
755  printf(" Rmin2 = %11.5f\n", fRmin2);
756  printf(" Rmax2 = %11.5f\n", fRmax2);
757  printf(" Bounding box:\n");
759 }
760 
761 ////////////////////////////////////////////////////////////////////////////////
762 /// Creates a TBuffer3D describing *this* shape.
763 /// Coordinates are in local reference frame.
764 
766 {
768  Int_t nbPnts = 4*n;
769  Int_t nbSegs = 8*n;
770  Int_t nbPols = 4*n;
772  nbPnts, 3*nbPnts,
773  nbSegs, 3*nbSegs,
774  nbPols, 6*nbPols);
775 
776  if (buff)
777  {
778  SetPoints(buff->fPnts);
779  SetSegsAndPols(*buff);
780  }
781  return buff;
782 }
783 
784 ////////////////////////////////////////////////////////////////////////////////
785 /// Fill TBuffer3D structure for segments and polygons.
786 
788 {
789  Int_t i,j;
791  Int_t c = GetBasicColor();
792 
793  for (i = 0; i < 4; i++) {
794  for (j = 0; j < n; j++) {
795  buffer.fSegs[(i*n+j)*3 ] = c;
796  buffer.fSegs[(i*n+j)*3+1] = i*n+j;
797  buffer.fSegs[(i*n+j)*3+2] = i*n+j+1;
798  }
799  buffer.fSegs[(i*n+j-1)*3+2] = i*n;
800  }
801  for (i = 4; i < 6; i++) {
802  for (j = 0; j < n; j++) {
803  buffer.fSegs[(i*n+j)*3 ] = c+1;
804  buffer.fSegs[(i*n+j)*3+1] = (i-4)*n+j;
805  buffer.fSegs[(i*n+j)*3+2] = (i-2)*n+j;
806  }
807  }
808  for (i = 6; i < 8; i++) {
809  for (j = 0; j < n; j++) {
810  buffer.fSegs[(i*n+j)*3 ] = c;
811  buffer.fSegs[(i*n+j)*3+1] = 2*(i-6)*n+j;
812  buffer.fSegs[(i*n+j)*3+2] = (2*(i-6)+1)*n+j;
813  }
814  }
815 
816  Int_t indx = 0;
817  i=0;
818  for (j = 0; j < n; j++) {
819  indx = 6*(i*n+j);
820  buffer.fPols[indx ] = c;
821  buffer.fPols[indx+1] = 4;
822  buffer.fPols[indx+5] = i*n+j;
823  buffer.fPols[indx+4] = (4+i)*n+j;
824  buffer.fPols[indx+3] = (2+i)*n+j;
825  buffer.fPols[indx+2] = (4+i)*n+j+1;
826  }
827  buffer.fPols[indx+2] = (4+i)*n;
828  i=1;
829  for (j = 0; j < n; j++) {
830  indx = 6*(i*n+j);
831  buffer.fPols[indx ] = c;
832  buffer.fPols[indx+1] = 4;
833  buffer.fPols[indx+2] = i*n+j;
834  buffer.fPols[indx+3] = (4+i)*n+j;
835  buffer.fPols[indx+4] = (2+i)*n+j;
836  buffer.fPols[indx+5] = (4+i)*n+j+1;
837  }
838  buffer.fPols[indx+5] = (4+i)*n;
839  i=2;
840  for (j = 0; j < n; j++) {
841  indx = 6*(i*n+j);
842  buffer.fPols[indx ] = c+i;
843  buffer.fPols[indx+1] = 4;
844  buffer.fPols[indx+2] = (i-2)*2*n+j;
845  buffer.fPols[indx+3] = (4+i)*n+j;
846  buffer.fPols[indx+4] = ((i-2)*2+1)*n+j;
847  buffer.fPols[indx+5] = (4+i)*n+j+1;
848  }
849  buffer.fPols[indx+5] = (4+i)*n;
850  i=3;
851  for (j = 0; j < n; j++) {
852  indx = 6*(i*n+j);
853  buffer.fPols[indx ] = c+i;
854  buffer.fPols[indx+1] = 4;
855  buffer.fPols[indx+5] = (i-2)*2*n+j;
856  buffer.fPols[indx+4] = (4+i)*n+j;
857  buffer.fPols[indx+3] = ((i-2)*2+1)*n+j;
858  buffer.fPols[indx+2] = (4+i)*n+j+1;
859  }
860  buffer.fPols[indx+2] = (4+i)*n;
861 }
862 
863 ////////////////////////////////////////////////////////////////////////////////
864 /// computes the closest distance from given point to this shape, according
865 /// to option. The matching point on the shape is stored in spoint.
866 
867 Double_t TGeoCone::Safety(const Double_t *point, Bool_t in) const
868 {
869  Double_t saf[4];
870  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
871  saf[0] = TGeoShape::SafetySeg(r,point[2], fRmin1, -fDz, fRmax1, -fDz, !in);
872  saf[1] = TGeoShape::SafetySeg(r,point[2], fRmax2, fDz, fRmin2, fDz, !in);
873  saf[2] = TGeoShape::SafetySeg(r,point[2], fRmin2, fDz, fRmin1, -fDz, !in);
874  saf[3] = TGeoShape::SafetySeg(r,point[2], fRmax1, -fDz, fRmax2, fDz, !in);
875  Double_t safety = saf[TMath::LocMin(4,saf)];
876  if (safety>1.E20) safety = 0.;
877  return safety;
878 }
879 
880 ////////////////////////////////////////////////////////////////////////////////
881 /// computes the closest distance from given point to this shape, according
882 /// to option. The matching point on the shape is stored in spoint.
883 
885  Double_t rmin2, Double_t rmax2, Int_t skipz)
886 {
887  Double_t saf[4];
888  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
889 // Double_t rin = tg1*point[2]+ro1;
890 // Double_t rout = tg2*point[2]+ro2;
891  switch (skipz) {
892  case 1: // skip lower Z plane
893  saf[0] = TGeoShape::Big();
894  saf[1] = TGeoShape::SafetySeg(r,point[2], rmax2, dz, rmin2, dz, !in);
895  break;
896  case 2: // skip upper Z plane
897  saf[0] = TGeoShape::SafetySeg(r,point[2], rmin1, -dz, rmax1, -dz, !in);
898  saf[1] = TGeoShape::Big();
899  break;
900  case 3: // skip both
901  saf[0] = saf[1] = TGeoShape::Big();
902  break;
903  default:
904  saf[0] = TGeoShape::SafetySeg(r,point[2], rmin1, -dz, rmax1, -dz, !in);
905  saf[1] = TGeoShape::SafetySeg(r,point[2], rmax2, dz, rmin2, dz, !in);
906  }
907  // Safety to inner part
908  if (rmin1>0 || rmin2>0)
909  saf[2] = TGeoShape::SafetySeg(r,point[2], rmin2, dz, rmin1, -dz, !in);
910  else
911  saf[2] = TGeoShape::Big();
912  saf[3] = TGeoShape::SafetySeg(r,point[2], rmax1, -dz, rmax2, dz, !in);
913  return saf[TMath::LocMin(4,saf)];
914 }
915 
916 ////////////////////////////////////////////////////////////////////////////////
917 /// Save a primitive as a C++ statement(s) on output stream "out".
918 
919 void TGeoCone::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
920 {
921  if (TObject::TestBit(kGeoSavePrimitive)) return;
922  out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
923  out << " dz = " << fDz << ";" << std::endl;
924  out << " rmin1 = " << fRmin1 << ";" << std::endl;
925  out << " rmax1 = " << fRmax1 << ";" << std::endl;
926  out << " rmin2 = " << fRmin2 << ";" << std::endl;
927  out << " rmax2 = " << fRmax2 << ";" << std::endl;
928  out << " TGeoShape *" << GetPointerName() << " = new TGeoCone(\"" << GetName() << "\", dz,rmin1,rmax1,rmin2,rmax2);" << std::endl;
930 }
931 
932 ////////////////////////////////////////////////////////////////////////////////
933 /// Set cone dimensions.
934 
936  Double_t rmin2, Double_t rmax2)
937 {
938  if (rmin1>=0) {
939  if (rmax1>0) {
940  if (rmin1<=rmax1) {
941  // normal rmin/rmax
942  fRmin1 = rmin1;
943  fRmax1 = rmax1;
944  } else {
945  fRmin1 = rmax1;
946  fRmax1 = rmin1;
947  Warning("SetConeDimensions", "rmin1>rmax1 Switch rmin1<->rmax1");
949  }
950  } else {
951  // run-time
952  fRmin1 = rmin1;
953  fRmax1 = rmax1;
954  }
955  } else {
956  // run-time
957  fRmin1 = rmin1;
958  fRmax1 = rmax1;
959  }
960  if (rmin2>=0) {
961  if (rmax2>0) {
962  if (rmin2<=rmax2) {
963  // normal rmin/rmax
964  fRmin2 = rmin2;
965  fRmax2 = rmax2;
966  } else {
967  fRmin2 = rmax2;
968  fRmax2 = rmin2;
969  Warning("SetConeDimensions", "rmin2>rmax2 Switch rmin2<->rmax2");
971  }
972  } else {
973  // run-time
974  fRmin2 = rmin2;
975  fRmax2 = rmax2;
976  }
977  } else {
978  // run-time
979  fRmin2 = rmin2;
980  fRmax2 = rmax2;
981  }
982 
983  fDz = dz;
984 }
985 
986 ////////////////////////////////////////////////////////////////////////////////
987 /// Set cone dimensions from an array.
988 
990 {
991  Double_t dz = param[0];
992  Double_t rmin1 = param[1];
993  Double_t rmax1 = param[2];
994  Double_t rmin2 = param[3];
995  Double_t rmax2 = param[4];
996  SetConeDimensions(dz, rmin1, rmax1, rmin2, rmax2);
997 }
998 
999 ////////////////////////////////////////////////////////////////////////////////
1000 /// Create cone mesh points.
1001 
1003 {
1004  Double_t dz, phi, dphi;
1005  Int_t j, n;
1006 
1007  n = gGeoManager->GetNsegments();
1008  dphi = 360./n;
1009  dz = fDz;
1010  Int_t indx = 0;
1011 
1012  if (points) {
1013  for (j = 0; j < n; j++) {
1014  phi = j*dphi*TMath::DegToRad();
1015  points[indx++] = fRmin1 * TMath::Cos(phi);
1016  points[indx++] = fRmin1 * TMath::Sin(phi);
1017  points[indx++] = -dz;
1018  }
1019 
1020  for (j = 0; j < n; j++) {
1021  phi = j*dphi*TMath::DegToRad();
1022  points[indx++] = fRmax1 * TMath::Cos(phi);
1023  points[indx++] = fRmax1 * TMath::Sin(phi);
1024  points[indx++] = -dz;
1025  }
1026 
1027  for (j = 0; j < n; j++) {
1028  phi = j*dphi*TMath::DegToRad();
1029  points[indx++] = fRmin2 * TMath::Cos(phi);
1030  points[indx++] = fRmin2 * TMath::Sin(phi);
1031  points[indx++] = dz;
1032  }
1033 
1034  for (j = 0; j < n; j++) {
1035  phi = j*dphi*TMath::DegToRad();
1036  points[indx++] = fRmax2 * TMath::Cos(phi);
1037  points[indx++] = fRmax2 * TMath::Sin(phi);
1038  points[indx++] = dz;
1039  }
1040  }
1041 }
1042 
1043 ////////////////////////////////////////////////////////////////////////////////
1044 /// Create cone mesh points.
1045 
1047 {
1048  Double_t dz, phi, dphi;
1049  Int_t j, n;
1050 
1051  n = gGeoManager->GetNsegments();
1052  dphi = 360./n;
1053  dz = fDz;
1054  Int_t indx = 0;
1055 
1056  if (points) {
1057  for (j = 0; j < n; j++) {
1058  phi = j*dphi*TMath::DegToRad();
1059  points[indx++] = fRmin1 * TMath::Cos(phi);
1060  points[indx++] = fRmin1 * TMath::Sin(phi);
1061  points[indx++] = -dz;
1062  }
1063 
1064  for (j = 0; j < n; j++) {
1065  phi = j*dphi*TMath::DegToRad();
1066  points[indx++] = fRmax1 * TMath::Cos(phi);
1067  points[indx++] = fRmax1 * TMath::Sin(phi);
1068  points[indx++] = -dz;
1069  }
1070 
1071  for (j = 0; j < n; j++) {
1072  phi = j*dphi*TMath::DegToRad();
1073  points[indx++] = fRmin2 * TMath::Cos(phi);
1074  points[indx++] = fRmin2 * TMath::Sin(phi);
1075  points[indx++] = dz;
1076  }
1077 
1078  for (j = 0; j < n; j++) {
1079  phi = j*dphi*TMath::DegToRad();
1080  points[indx++] = fRmax2 * TMath::Cos(phi);
1081  points[indx++] = fRmax2 * TMath::Sin(phi);
1082  points[indx++] = dz;
1083  }
1084  }
1085 }
1086 
1087 ////////////////////////////////////////////////////////////////////////////////
1088 /// Returns numbers of vertices, segments and polygons composing the shape mesh.
1089 
1090 void TGeoCone::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
1091 {
1093  nvert = n*4;
1094  nsegs = n*8;
1095  npols = n*4;
1096 }
1097 
1098 ////////////////////////////////////////////////////////////////////////////////
1099 /// Return number of vertices of the mesh representation
1100 
1102 {
1104  Int_t numPoints = n*4;
1105  return numPoints;
1106 }
1107 
1108 ////////////////////////////////////////////////////////////////////////////////
1109 /// Fill size of this 3-D object
1110 
1112 {
1113 }
1114 
1115 ////////////////////////////////////////////////////////////////////////////////
1116 /// Fills a static 3D buffer and returns a reference.
1117 
1118 const TBuffer3D & TGeoCone::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
1119 {
1120  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
1121 
1122  TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
1123 
1124  if (reqSections & TBuffer3D::kRawSizes) {
1126  Int_t nbPnts = 4*n;
1127  Int_t nbSegs = 8*n;
1128  Int_t nbPols = 4*n;
1129  if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
1131  }
1132  }
1133 
1134  // TODO: Can we push this as common down to TGeoShape?
1135  if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
1136  SetPoints(buffer.fPnts);
1137  if (!buffer.fLocalFrame) {
1138  TransformPoints(buffer.fPnts, buffer.NbPnts());
1139  }
1140 
1141  SetSegsAndPols(buffer);
1143  }
1144 
1145  return buffer;
1146 }
1147 
1148 ////////////////////////////////////////////////////////////////////////////////
1149 /// Check the inside status for each of the points in the array.
1150 /// Input: Array of point coordinates + vector size
1151 /// Output: Array of Booleans for the inside of each point
1152 
1153 void TGeoCone::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
1154 {
1155  for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
1156 }
1157 
1158 ////////////////////////////////////////////////////////////////////////////////
1159 /// Compute the normal for an array o points so that norm.dot.dir is positive
1160 /// Input: Arrays of point coordinates and directions + vector size
1161 /// Output: Array of normal directions
1162 
1163 void TGeoCone::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
1164 {
1165  for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
1166 }
1167 
1168 ////////////////////////////////////////////////////////////////////////////////
1169 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
1170 
1171 void TGeoCone::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1172 {
1173  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
1174 }
1175 
1176 ////////////////////////////////////////////////////////////////////////////////
1177 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
1178 
1179 void TGeoCone::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1180 {
1181  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
1182 }
1183 
1184 ////////////////////////////////////////////////////////////////////////////////
1185 /// Compute safe distance from each of the points in the input array.
1186 /// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1187 /// Output: Safety values
1188 
1189 void TGeoCone::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
1190 {
1191  for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
1192 }
1193 
1195 
1196 ////////////////////////////////////////////////////////////////////////////////
1197 /// Default constructor
1198 
1200  :TGeoCone(),
1201  fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1202 {
1204  fPhi1 = fPhi2 = 0.0;
1205 }
1206 
1207 ////////////////////////////////////////////////////////////////////////////////
1208 /// Default constructor specifying minimum and maximum radius
1209 
1211  Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2)
1212  :TGeoCone(dz, rmin1, rmax1, rmin2, rmax2),
1213  fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1214 
1215 {
1217  SetConsDimensions(dz, rmin1, rmax1, rmin2, rmax2, phi1, phi2);
1218  ComputeBBox();
1219 }
1220 
1221 ////////////////////////////////////////////////////////////////////////////////
1222 /// Default constructor specifying minimum and maximum radius
1223 
1225  Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2)
1226  :TGeoCone(name, dz, rmin1, rmax1, rmin2, rmax2),
1227  fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1228 {
1230  SetConsDimensions(dz, rmin1, rmax1, rmin2, rmax2, phi1, phi2);
1231  ComputeBBox();
1232 }
1233 
1234 ////////////////////////////////////////////////////////////////////////////////
1235 /// Default constructor specifying minimum and maximum radius
1236 /// - param[0] = dz
1237 /// - param[1] = Rmin1
1238 /// - param[2] = Rmax1
1239 /// - param[3] = Rmin2
1240 /// - param[4] = Rmax2
1241 /// - param[5] = phi1
1242 /// - param[6] = phi2
1243 
1245  :TGeoCone(0,0,0,0,0),
1246  fPhi1(0.), fPhi2(0.), fS1(0.), fC1(0.), fS2(0.), fC2(0.), fSm(0.), fCm(0.), fCdfi(0.)
1247 {
1249  SetDimensions(param);
1250  ComputeBBox();
1251 }
1252 
1253 ////////////////////////////////////////////////////////////////////////////////
1254 /// destructor
1255 
1257 {
1258 }
1259 
1260 ////////////////////////////////////////////////////////////////////////////////
1261 /// Function called after streaming an object of this class.
1262 
1264 {
1265  InitTrigonometry();
1266 }
1267 
1268 ////////////////////////////////////////////////////////////////////////////////
1269 /// Init frequently used trigonometric values
1270 
1272 {
1273  Double_t phi1 = fPhi1*TMath::DegToRad();
1274  Double_t phi2 = fPhi2*TMath::DegToRad();
1275  fC1 = TMath::Cos(phi1);
1276  fS1 = TMath::Sin(phi1);
1277  fC2 = TMath::Cos(phi2);
1278  fS2 = TMath::Sin(phi2);
1279  Double_t fio = 0.5*(phi1+phi2);
1280  fCm = TMath::Cos(fio);
1281  fSm = TMath::Sin(fio);
1282  Double_t dfi = 0.5*(phi2-phi1);
1283  fCdfi = TMath::Cos(dfi);
1284 }
1285 
1286 ////////////////////////////////////////////////////////////////////////////////
1287 /// Computes capacity of the shape in [length^3]
1288 
1290 {
1292 }
1293 
1294 ////////////////////////////////////////////////////////////////////////////////
1295 /// Computes capacity of the shape in [length^3]
1296 
1298 {
1299  Double_t capacity = (TMath::Abs(phi2-phi1)*TMath::DegToRad()*dz/3.)*
1300  (rmax1*rmax1+rmax2*rmax2+rmax1*rmax2-
1301  rmin1*rmin1-rmin2*rmin2-rmin1*rmin2);
1302  return capacity;
1303 }
1304 
1305 ////////////////////////////////////////////////////////////////////////////////
1306 /// compute bounding box of the tube segment
1307 
1309 {
1310  Double_t rmin, rmax;
1311  rmin = TMath::Min(fRmin1, fRmin2);
1312  rmax = TMath::Max(fRmax1, fRmax2);
1313 
1314  Double_t xc[4];
1315  Double_t yc[4];
1316  xc[0] = rmax*fC1;
1317  yc[0] = rmax*fS1;
1318  xc[1] = rmax*fC2;
1319  yc[1] = rmax*fS2;
1320  xc[2] = rmin*fC1;
1321  yc[2] = rmin*fS1;
1322  xc[3] = rmin*fC2;
1323  yc[3] = rmin*fS2;
1324 
1325  Double_t xmin = xc[TMath::LocMin(4, &xc[0])];
1326  Double_t xmax = xc[TMath::LocMax(4, &xc[0])];
1327  Double_t ymin = yc[TMath::LocMin(4, &yc[0])];
1328  Double_t ymax = yc[TMath::LocMax(4, &yc[0])];
1329 
1330  Double_t dp = fPhi2-fPhi1;
1331  Double_t ddp = -fPhi1;
1332  if (ddp<0) ddp+= 360;
1333  if (ddp<=dp) xmax = rmax;
1334  ddp = 90-fPhi1;
1335  if (ddp<0) ddp+= 360;
1336  if (ddp<=dp) ymax = rmax;
1337  ddp = 180-fPhi1;
1338  if (ddp<0) ddp+= 360;
1339  if (ddp<=dp) xmin = -rmax;
1340  ddp = 270-fPhi1;
1341  if (ddp<0) ddp+= 360;
1342  if (ddp<=dp) ymin = -rmax;
1343  fOrigin[0] = (xmax+xmin)/2;
1344  fOrigin[1] = (ymax+ymin)/2;
1345  fOrigin[2] = 0;
1346  fDX = (xmax-xmin)/2;
1347  fDY = (ymax-ymin)/2;
1348  fDZ = fDz;
1349 }
1350 
1351 ////////////////////////////////////////////////////////////////////////////////
1352 /// Compute normal to closest surface from POINT.
1353 
1354 void TGeoConeSeg::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
1355 {
1356  Double_t saf[3];
1357  Double_t ro1 = 0.5*(fRmin1+fRmin2);
1358  Double_t tg1 = 0.5*(fRmin2-fRmin1)/fDz;
1359  Double_t cr1 = 1./TMath::Sqrt(1.+tg1*tg1);
1360  Double_t ro2 = 0.5*(fRmax1+fRmax2);
1361  Double_t tg2 = 0.5*(fRmax2-fRmax1)/fDz;
1362  Double_t cr2 = 1./TMath::Sqrt(1.+tg2*tg2);
1363 
1364  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
1365  Double_t rin = tg1*point[2]+ro1;
1366  Double_t rout = tg2*point[2]+ro2;
1367  saf[0] = TMath::Abs(fDz-TMath::Abs(point[2]));
1368  saf[1] = (ro1>0)?(TMath::Abs((r-rin)*cr1)):TGeoShape::Big();
1369  saf[2] = TMath::Abs((rout-r)*cr2);
1370  Int_t i = TMath::LocMin(3,saf);
1371  if (((fPhi2-fPhi1)<360.) && TGeoShape::IsCloseToPhi(saf[i], point,fC1,fS1,fC2,fS2)) {
1372  TGeoShape::NormalPhi(point,dir,norm,fC1,fS1,fC2,fS2);
1373  return;
1374  }
1375  if (i==0) {
1376  norm[0] = norm[1] = 0.;
1377  norm[2] = TMath::Sign(1.,dir[2]);
1378  return;
1379  }
1380 
1381  Double_t phi = TMath::ATan2(point[1],point[0]);
1382  Double_t cphi = TMath::Cos(phi);
1383  Double_t sphi = TMath::Sin(phi);
1384 
1385  if (i==1) {
1386  norm[0] = cr1*cphi;
1387  norm[1] = cr1*sphi;
1388  norm[2] = -tg1*cr1;
1389  } else {
1390  norm[0] = cr2*cphi;
1391  norm[1] = cr2*sphi;
1392  norm[2] = -tg2*cr2;
1393  }
1394 
1395  if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
1396  norm[0] = -norm[0];
1397  norm[1] = -norm[1];
1398  norm[2] = -norm[2];
1399  }
1400 }
1401 
1402 ////////////////////////////////////////////////////////////////////////////////
1403 /// Compute normal to closest surface from POINT.
1404 
1405 void TGeoConeSeg::ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm,
1406  Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2,
1408 {
1409  Double_t saf[2];
1410  Double_t ro1 = 0.5*(rmin1+rmin2);
1411  Double_t tg1 = 0.5*(rmin2-rmin1)/dz;
1412  Double_t cr1 = 1./TMath::Sqrt(1.+tg1*tg1);
1413  Double_t ro2 = 0.5*(rmax1+rmax2);
1414  Double_t tg2 = 0.5*(rmax2-rmax1)/dz;
1415  Double_t cr2 = 1./TMath::Sqrt(1.+tg2*tg2);
1416 
1417  Double_t r=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
1418  Double_t rin = tg1*point[2]+ro1;
1419  Double_t rout = tg2*point[2]+ro2;
1420  saf[0] = (ro1>0)?(TMath::Abs((r-rin)*cr1)):TGeoShape::Big();
1421  saf[1] = TMath::Abs((rout-r)*cr2);
1422  Int_t i = TMath::LocMin(2,saf);
1423  if (TGeoShape::IsCloseToPhi(saf[i], point,c1,s1,c2,s2)) {
1424  TGeoShape::NormalPhi(point,dir,norm,c1,s1,c2,s2);
1425  return;
1426  }
1427 
1428  Double_t phi = TMath::ATan2(point[1],point[0]);
1429  Double_t cphi = TMath::Cos(phi);
1430  Double_t sphi = TMath::Sin(phi);
1431 
1432  if (i==0) {
1433  norm[0] = cr1*cphi;
1434  norm[1] = cr1*sphi;
1435  norm[2] = -tg1*cr1;
1436  } else {
1437  norm[0] = cr2*cphi;
1438  norm[1] = cr2*sphi;
1439  norm[2] = -tg2*cr2;
1440  }
1441 
1442  if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
1443  norm[0] = -norm[0];
1444  norm[1] = -norm[1];
1445  norm[2] = -norm[2];
1446  }
1447 }
1448 
1449 ////////////////////////////////////////////////////////////////////////////////
1450 /// test if point is inside this sphere
1451 
1453 {
1454  if (!TGeoCone::Contains(point)) return kFALSE;
1455  Double_t dphi = fPhi2 - fPhi1;
1456  if (dphi >= 360.) return kTRUE;
1457  Double_t phi = TMath::ATan2(point[1], point[0]) * TMath::RadToDeg();
1458  if (phi < 0 ) phi+=360.;
1459  Double_t ddp = phi-fPhi1;
1460  if (ddp < 0) ddp+=360.;
1461 // if (ddp > 360) ddp-=360;
1462  if (ddp > dphi) return kFALSE;
1463  return kTRUE;
1464 }
1465 
1466 ////////////////////////////////////////////////////////////////////////////////
1467 /// Static method to compute distance to a conical surface with :
1468 /// - r1, z1 - radius and Z position of lower base
1469 /// - r2, z2 - radius and Z position of upper base
1470 /// - phi1, phi2 - phi limits
1471 
1473 {
1474  Double_t dz = z2-z1;
1475  if (dz<=0) {
1476  return TGeoShape::Big();
1477  }
1478 
1479  Double_t dphi = phi2 - phi1;
1480  Bool_t hasphi = kTRUE;
1481  if (dphi >= 360.) hasphi=kFALSE;
1482  if (dphi < 0) dphi+=360.;
1483 // printf("phi1=%f phi2=%f dphi=%f\n", phi1, phi2, dphi);
1484 
1485  Double_t ro0 = 0.5*(r1+r2);
1486  Double_t fz = (r2-r1)/dz;
1487  Double_t r0sq = point[0]*point[0] + point[1]*point[1];
1488  Double_t rc = ro0 + fz*(point[2]-0.5*(z1+z2));
1489 
1490  Double_t a = dir[0]*dir[0] + dir[1]*dir[1] - fz*fz*dir[2]*dir[2];
1491  Double_t b = point[0]*dir[0] + point[1]*dir[1] - fz*rc*dir[2];
1492  Double_t c = r0sq - rc*rc;
1493 
1494  if (a==0) return TGeoShape::Big();
1495  a = 1./a;
1496  b *= a;
1497  c *= a;
1498  Double_t delta = b*b - c;
1499  if (delta<0) return TGeoShape::Big();
1500  delta = TMath::Sqrt(delta);
1501 
1502  Double_t snxt = -b-delta;
1503  Double_t ptnew[3];
1504  Double_t ddp, phi;
1505  if (snxt>0) {
1506  // check Z range
1507  ptnew[2] = point[2] + snxt*dir[2];
1508  if (((ptnew[2]-z1)*(ptnew[2]-z2)) < 0) {
1509  // check phi range
1510  if (!hasphi) return snxt;
1511  ptnew[0] = point[0] + snxt*dir[0];
1512  ptnew[1] = point[1] + snxt*dir[1];
1513  phi = TMath::ATan2(ptnew[1], ptnew[0]) * TMath::RadToDeg();
1514  if (phi < 0 ) phi+=360.;
1515  ddp = phi-phi1;
1516  if (ddp < 0) ddp+=360.;
1517  // printf("snxt1=%f phi=%f ddp=%f\n", snxt, phi, ddp);
1518  if (ddp<=dphi) return snxt;
1519  }
1520  }
1521  snxt = -b+delta;
1522  if (snxt>0) {
1523  // check Z range
1524  ptnew[2] = point[2] + snxt*dir[2];
1525  if (((ptnew[2]-z1)*(ptnew[2]-z2)) < 0) {
1526  // check phi range
1527  if (!hasphi) return snxt;
1528  ptnew[0] = point[0] + snxt*dir[0];
1529  ptnew[1] = point[1] + snxt*dir[1];
1530  phi = TMath::ATan2(ptnew[1], ptnew[0]) * TMath::RadToDeg();
1531  if (phi < 0 ) phi+=360.;
1532  ddp = phi-phi1;
1533  if (ddp < 0) ddp+=360.;
1534  // printf("snxt2=%f phi=%f ddp=%f\n", snxt, phi, ddp);
1535  if (ddp<=dphi) return snxt;
1536  }
1537  }
1538  return TGeoShape::Big();
1539 }
1540 
1541 ////////////////////////////////////////////////////////////////////////////////
1542 /// compute distance from inside point to surface of the tube segment
1543 
1545  Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2,
1547 {
1548  if (dz<=0) return TGeoShape::Big();
1549  // Do Z
1550  Double_t scone = TGeoCone::DistFromInsideS(point,dir,dz,rmin1,rmax1,rmin2,rmax2);
1551  if (scone<=0) return 0.0;
1552  Double_t sfmin = TGeoShape::Big();
1553  Double_t rsq = point[0]*point[0]+point[1]*point[1];
1554  Double_t r = TMath::Sqrt(rsq);
1555  Double_t cpsi=point[0]*cm+point[1]*sm;
1556  if (cpsi>r*cdfi+TGeoShape::Tolerance()) {
1557  sfmin = TGeoShape::DistToPhiMin(point, dir, s1, c1, s2, c2, sm, cm);
1558  return TMath::Min(scone,sfmin);
1559  }
1560  // Point on the phi boundary or outside
1561  // which one: phi1 or phi2
1562  Double_t ddotn, xi, yi;
1563  if (TMath::Abs(point[1]-s1*r) < TMath::Abs(point[1]-s2*r)) {
1564  ddotn = s1*dir[0]-c1*dir[1];
1565  if (ddotn>=0) return 0.0;
1566  ddotn = -s2*dir[0]+c2*dir[1];
1567  if (ddotn<=0) return scone;
1568  sfmin = s2*point[0]-c2*point[1];
1569  if (sfmin<=0) return scone;
1570  sfmin /= ddotn;
1571  if (sfmin >= scone) return scone;
1572  xi = point[0]+sfmin*dir[0];
1573  yi = point[1]+sfmin*dir[1];
1574  if (yi*cm-xi*sm<0) return scone;
1575  return sfmin;
1576  }
1577  ddotn = -s2*dir[0]+c2*dir[1];
1578  if (ddotn>=0) return 0.0;
1579  ddotn = s1*dir[0]-c1*dir[1];
1580  if (ddotn<=0) return scone;
1581  sfmin = -s1*point[0]+c1*point[1];
1582  if (sfmin<=0) return scone;
1583  sfmin /= ddotn;
1584  if (sfmin >= scone) return scone;
1585  xi = point[0]+sfmin*dir[0];
1586  yi = point[1]+sfmin*dir[1];
1587  if (yi*cm-xi*sm>0) return scone;
1588  return sfmin;
1589 }
1590 
1591 ////////////////////////////////////////////////////////////////////////////////
1592 /// compute distance from inside point to surface of the tube segment
1593 
1594 Double_t TGeoConeSeg::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1595 {
1596  if (iact<3 && safe) {
1598  if (iact==0) return TGeoShape::Big();
1599  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
1600  }
1601  if ((fPhi2-fPhi1)>=360.) return TGeoCone::DistFromInsideS(point,dir,fDz,fRmin1,fRmax1,fRmin2,fRmax2);
1602 
1603  // compute distance to surface
1605 }
1606 
1607 ////////////////////////////////////////////////////////////////////////////////
1608 /// compute distance from outside point to surface of arbitrary tube
1609 
1611  Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2,
1613 {
1614  if (dz<=0) return TGeoShape::Big();
1615  Double_t r2, cpsi;
1616  // check Z planes
1617  Double_t xi, yi, zi;
1618  Double_t b,delta;
1619  zi = dz - TMath::Abs(point[2]);
1620  Double_t rin,rout;
1622  Double_t snxt=TGeoShape::Big();
1623  Bool_t in = kFALSE;
1624  Bool_t inz = (zi<0)?kFALSE:kTRUE;
1625  if (!inz) {
1626  if (point[2]*dir[2]>=0) return TGeoShape::Big();
1627  s = -zi/TMath::Abs(dir[2]);
1628  xi = point[0]+s*dir[0];
1629  yi = point[1]+s*dir[1];
1630  r2=xi*xi+yi*yi;
1631  if (dir[2]>0) {
1632  rin = rmin1;
1633  rout = rmax1;
1634  } else {
1635  rin = rmin2;
1636  rout = rmax2;
1637  }
1638  if ((rin*rin<=r2) && (r2<=rout*rout)) {
1639  cpsi=xi*cm+yi*sm;
1640  if (cpsi>=(cdfi*TMath::Sqrt(r2))) return s;
1641  }
1642  }
1643  Double_t zinv = 1./dz;
1644  Double_t rsq = point[0]*point[0]+point[1]*point[1];
1645  Double_t r = TMath::Sqrt(rsq);
1646  Double_t ro1=0.5*(rmin1+rmin2);
1647  Bool_t hasrmin = (ro1>0)?kTRUE:kFALSE;
1648  Double_t tg1 = 0.0;
1649  Bool_t inrmin = kFALSE;
1650  rin = 0.0;
1651  if (hasrmin) {
1652  tg1=0.5*(rmin2-rmin1)*zinv;
1653  rin = ro1+tg1*point[2];
1654  if (rsq > rin*(rin-TGeoShape::Tolerance())) inrmin=kTRUE;
1655  } else {
1656  inrmin = kTRUE;
1657  }
1658  Double_t ro2=0.5*(rmax1+rmax2);
1659  Double_t tg2=0.5*(rmax2-rmax1)*zinv;
1660  rout = ro2+tg2*point[2];
1661  Bool_t inrmax = kFALSE;
1662  if (r < rout+TGeoShape::Tolerance()) inrmax = kTRUE;
1663  Bool_t inphi = kFALSE;
1664  cpsi=point[0]*cm+point[1]*sm;
1665  if (cpsi>r*cdfi-TGeoShape::Tolerance()) inphi = kTRUE;
1666  in = inz & inrmin & inrmax & inphi;
1667  // If inside, we are most likely on a boundary within machine precision.
1668  if (in) {
1669  Double_t safphi = (cpsi-r*cdfi)*TMath::Sqrt(1.-cdfi*cdfi);
1670  Double_t safrmin = (hasrmin)?TMath::Abs(r-rin):(TGeoShape::Big());
1671  Double_t safrmax = TMath::Abs(r-rout);
1672  // check if on Z boundaries
1673  if (zi<safrmax && zi<safrmin && zi<safphi) {
1674  if (point[2]*dir[2]<0) return 0.0;
1675  return TGeoShape::Big();
1676  }
1677  // check if on Rmax boundary
1678  if (safrmax<safrmin && safrmax<safphi) {
1679  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]-tg2*dir[2]*r;
1680  if (ddotn<=0) return 0.0;
1681  return TGeoShape::Big();
1682  }
1683  // check if on phi boundary
1684  if (safphi<safrmin) {
1685  // We may cross again a phi of rmin boundary
1686  // check first if we are on phi1 or phi2
1687  Double_t un;
1688  if (point[0]*c1 + point[1]*s1 > point[0]*c2 + point[1]*s2) {
1689  un = dir[0]*s1-dir[1]*c1;
1690  if (un < 0) return 0.0;
1691  if (cdfi>=0) return TGeoShape::Big();
1692  un = -dir[0]*s2+dir[1]*c2;
1693  if (un<0) {
1694  s = -point[0]*s2+point[1]*c2;
1695  if (s>0) {
1696  s /= (-un);
1697  zi = point[2]+s*dir[2];
1698  if (TMath::Abs(zi)<=dz) {
1699  xi = point[0]+s*dir[0];
1700  yi = point[1]+s*dir[1];
1701  if ((yi*cm-xi*sm)>0) {
1702  r2=xi*xi+yi*yi;
1703  rin = ro1+tg1*zi;
1704  rout = ro2+tg2*zi;
1705  if ((rin*rin<=r2) && (rout*rout>=r2)) return s;
1706  }
1707  }
1708  }
1709  }
1710  } else {
1711  un = -dir[0]*s2+dir[1]*c2;
1712  if (un < 0) return 0.0;
1713  if (cdfi>=0) return TGeoShape::Big();
1714  un = dir[0]*s1-dir[1]*c1;
1715  if (un<0) {
1716  s = point[0]*s1-point[1]*c1;
1717  if (s>0) {
1718  s /= (-un);
1719  zi = point[2]+s*dir[2];
1720  if (TMath::Abs(zi)<=dz) {
1721  xi = point[0]+s*dir[0];
1722  yi = point[1]+s*dir[1];
1723  if ((yi*cm-xi*sm)<0) {
1724  r2=xi*xi+yi*yi;
1725  rin = ro1+tg1*zi;
1726  rout = ro2+tg2*zi;
1727  if ((rin*rin<=r2) && (rout*rout>=r2)) return s;
1728  }
1729  }
1730  }
1731  }
1732  }
1733  // We may also cross rmin, second solution coming from outside
1734  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]-tg1*dir[2]*r;
1735  if (ddotn>=0) return TGeoShape::Big();
1736  if (cdfi>=0) return TGeoShape::Big();
1737  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
1738  if (delta<0) return TGeoShape::Big();
1739  snxt = -b-delta;
1740  if (snxt<0) return TGeoShape::Big();
1741  snxt = -b+delta;
1742  zi = point[2]+snxt*dir[2];
1743  if (TMath::Abs(zi)>dz) return TGeoShape::Big();
1744  xi = point[0]+snxt*dir[0];
1745  yi = point[1]+snxt*dir[1];
1746  r2=xi*xi+yi*yi;
1747  cpsi=xi*cm+yi*sm;
1748  if (cpsi>=(cdfi*TMath::Sqrt(r2))) return snxt;
1749  return TGeoShape::Big();
1750  }
1751  // We are on rmin boundary: we may cross again rmin or a phi facette
1752  Double_t ddotn = point[0]*dir[0]+point[1]*dir[1]-tg1*dir[2]*r;
1753  if (ddotn>=0) return 0.0;
1754  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
1755  if (delta<0) return 0.0;
1756  snxt = -b+delta;
1757  if (snxt<0) return TGeoShape::Big();
1758  if (TMath::Abs(-b-delta)>snxt) return TGeoShape::Big();
1759  zi = point[2]+snxt*dir[2];
1760  if (TMath::Abs(zi)>dz) return TGeoShape::Big();
1761  // OK, we cross rmin at snxt - check if within phi range
1762  xi = point[0]+snxt*dir[0];
1763  yi = point[1]+snxt*dir[1];
1764  r2=xi*xi+yi*yi;
1765  cpsi=xi*cm+yi*sm;
1766  if (cpsi>=(cdfi*TMath::Sqrt(r2))) return snxt;
1767  // we cross rmin in the phi gap - we may cross a phi facette
1768  if (cdfi>=0) return TGeoShape::Big();
1769  Double_t un=-dir[0]*s1+dir[1]*c1;
1770  if (un > 0) {
1771  s=point[0]*s1-point[1]*c1;
1772  if (s>=0) {
1773  s /= un;
1774  zi=point[2]+s*dir[2];
1775  if (TMath::Abs(zi)<=dz) {
1776  xi=point[0]+s*dir[0];
1777  yi=point[1]+s*dir[1];
1778  if ((yi*cm-xi*sm)<=0) {
1779  r2=xi*xi+yi*yi;
1780  rin = ro1+tg1*zi;
1781  rout = ro2+tg2*zi;
1782  if ((rin*rin<=r2) && (rout*rout>=r2)) return s;
1783  }
1784  }
1785  }
1786  }
1787  un=dir[0]*s2-dir[1]*c2;
1788  if (un > 0) {
1789  s=(point[1]*c2-point[0]*s2)/un;
1790  if (s>=0) {
1791  zi=point[2]+s*dir[2];
1792  if (TMath::Abs(zi)<=dz) {
1793  xi=point[0]+s*dir[0];
1794  yi=point[1]+s*dir[1];
1795  if ((yi*cm-xi*sm)>=0) {
1796  r2=xi*xi+yi*yi;
1797  rin = ro1+tg1*zi;
1798  rout = ro2+tg2*zi;
1799  if ((rin*rin<=r2) && (rout*rout>=r2)) return s;
1800  }
1801  }
1802  }
1803  }
1804  return TGeoShape::Big();
1805  }
1806 
1807  // The point is really outside
1808  Double_t sr1 = TGeoShape::Big();
1809  if (!inrmax) {
1810  // check crossing with outer cone
1811  TGeoCone::DistToCone(point, dir, dz, rmax1, rmax2, b, delta);
1812  if (delta>=0) {
1813  s = -b-delta;
1814  if (s>0) {
1815  zi=point[2]+s*dir[2];
1816  if (TMath::Abs(zi)<=dz) {
1817  xi=point[0]+s*dir[0];
1818  yi=point[1]+s*dir[1];
1819  r2=xi*xi+yi*yi;
1820  cpsi=xi*cm+yi*sm;
1821  if (cpsi>=(cdfi*TMath::Sqrt(r2))) return s; // rmax crossing
1822  }
1823  }
1824  s = -b+delta;
1825  if (s>0) {
1826  zi=point[2]+s*dir[2];
1827  if (TMath::Abs(zi)<=dz) {
1828  xi=point[0]+s*dir[0];
1829  yi=point[1]+s*dir[1];
1830  r2=xi*xi+yi*yi;
1831  cpsi=xi*cm+yi*sm;
1832  if (cpsi>=(cdfi*TMath::Sqrt(r2))) sr1=s;
1833  }
1834  }
1835  }
1836  }
1837  // check crossing with inner cone
1838  Double_t sr2 = TGeoShape::Big();
1839  TGeoCone::DistToCone(point, dir, dz, rmin1, rmin2, b, delta);
1840  if (delta>=0) {
1841  s = -b-delta;
1842  if (s>0) {
1843  zi=point[2]+s*dir[2];
1844  if (TMath::Abs(zi)<=dz) {
1845  xi=point[0]+s*dir[0];
1846  yi=point[1]+s*dir[1];
1847  r2=xi*xi+yi*yi;
1848  cpsi=xi*cm+yi*sm;
1849  if (cpsi>=(cdfi*TMath::Sqrt(r2))) sr2=s;
1850  }
1851  }
1852  if (sr2>1E10) {
1853  s = -b+delta;
1854  if (s>0) {
1855  zi=point[2]+s*dir[2];
1856  if (TMath::Abs(zi)<=dz) {
1857  xi=point[0]+s*dir[0];
1858  yi=point[1]+s*dir[1];
1859  r2=xi*xi+yi*yi;
1860  cpsi=xi*cm+yi*sm;
1861  if (cpsi>=(cdfi*TMath::Sqrt(r2))) sr2=s;
1862  }
1863  }
1864  }
1865  }
1866  snxt = TMath::Min(sr1,sr2);
1867  // Check phi crossing
1868  s = TGeoShape::DistToPhiMin(point,dir,s1,c1,s2,c2,sm,cm,kFALSE);
1869  if (s>snxt) return snxt;
1870  zi=point[2]+s*dir[2];
1871  if (TMath::Abs(zi)>dz) return snxt;
1872  xi=point[0]+s*dir[0];
1873  yi=point[1]+s*dir[1];
1874  r2=xi*xi+yi*yi;
1875  rout = ro2+tg2*zi;
1876  if (r2>rout*rout) return snxt;
1877  rin = ro1+tg1*zi;
1878  if (r2>=rin*rin) return s; // phi crossing
1879  return snxt;
1880 }
1881 
1882 ////////////////////////////////////////////////////////////////////////////////
1883 /// compute distance from outside point to surface of the tube
1884 
1885 Double_t TGeoConeSeg::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1886 {
1887  // compute safe radius
1888  if (iact<3 && safe) {
1889  *safe = Safety(point, kFALSE);
1890  if (iact==0) return TGeoShape::Big();
1891  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
1892  }
1893  // Check if the bounding box is crossed within the requested distance
1894  Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
1895  if (sdist>=step) return TGeoShape::Big();
1896  if ((fPhi2-fPhi1)>=360.) return TGeoCone::DistFromOutsideS(point,dir,fDz,fRmin1,fRmax1,fRmin2,fRmax2);
1898 }
1899 
1900 ////////////////////////////////////////////////////////////////////////////////
1901 /// compute closest distance from point px,py to each corner
1902 
1904 {
1906  const Int_t numPoints = 4*n;
1907  return ShapeDistancetoPrimitive(numPoints, px, py);
1908 }
1909 
1910 ////////////////////////////////////////////////////////////////////////////////
1911 /// Divide this cone segment shape belonging to volume "voldiv" into ndiv volumes
1912 /// called divname, from start position with the given step. Returns pointer
1913 /// to created division cell volume in case of Z divisions. For Z division
1914 /// creates all volumes with different shapes and returns pointer to volume that
1915 /// was divided. In case a wrong division axis is supplied, returns pointer to
1916 /// volume that was divided.
1917 
1918 TGeoVolume *TGeoConeSeg::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
1919  Double_t start, Double_t step)
1920 {
1921  TGeoShape *shape; //--- shape to be created
1922  TGeoVolume *vol; //--- division volume to be created
1923  TGeoVolumeMulti *vmulti; //--- generic divided volume
1924  TGeoPatternFinder *finder; //--- finder to be attached
1925  TString opt = ""; //--- option to be attached
1926  Double_t dphi;
1927  Int_t id;
1928  Double_t end = start+ndiv*step;
1929  switch (iaxis) {
1930  case 1: //--- R division
1931  Error("Divide","division of a cone segment on R not implemented");
1932  return 0;
1933  case 2: //--- Phi division
1934  dphi = fPhi2-fPhi1;
1935  if (dphi<0) dphi+=360.;
1936  finder = new TGeoPatternCylPhi(voldiv, ndiv, start, end);
1937  voldiv->SetFinder(finder);
1938  finder->SetDivIndex(voldiv->GetNdaughters());
1939  shape = new TGeoConeSeg(fDz, fRmin1, fRmax1, fRmin2, fRmax2, -step/2, step/2);
1940  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1941  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1942  vmulti->AddVolume(vol);
1943  opt = "Phi";
1944  for (id=0; id<ndiv; id++) {
1945  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
1946  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1947  }
1948  return vmulti;
1949  case 3: //--- Z division
1950  finder = new TGeoPatternZ(voldiv, ndiv, start, end);
1951  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1952  voldiv->SetFinder(finder);
1953  finder->SetDivIndex(voldiv->GetNdaughters());
1954  for (id=0; id<ndiv; id++) {
1955  Double_t z1 = start+id*step;
1956  Double_t z2 = start+(id+1)*step;
1957  Double_t rmin1n = 0.5*(fRmin1*(fDz-z1)+fRmin2*(fDz+z1))/fDz;
1958  Double_t rmax1n = 0.5*(fRmax1*(fDz-z1)+fRmax2*(fDz+z1))/fDz;
1959  Double_t rmin2n = 0.5*(fRmin1*(fDz-z2)+fRmin2*(fDz+z2))/fDz;
1960  Double_t rmax2n = 0.5*(fRmax1*(fDz-z2)+fRmax2*(fDz+z2))/fDz;
1961  shape = new TGeoConeSeg(step/2, rmin1n, rmax1n, rmin2n, rmax2n, fPhi1, fPhi2);
1962  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1963  vmulti->AddVolume(vol);
1964  opt = "Z";
1965  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
1966  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1967  }
1968  return vmulti;
1969  default:
1970  Error("Divide", "Wrong axis type for division");
1971  return 0;
1972  }
1973 }
1974 
1975 ////////////////////////////////////////////////////////////////////////////////
1976 /// Get range of shape for a given axis.
1977 
1979 {
1980  xlo = 0;
1981  xhi = 0;
1982  Double_t dx = 0;
1983  switch (iaxis) {
1984  case 2:
1985  xlo = fPhi1;
1986  xhi = fPhi2;
1987  dx = xhi-xlo;
1988  return dx;
1989  case 3:
1990  xlo = -fDz;
1991  xhi = fDz;
1992  dx = xhi-xlo;
1993  return dx;
1994  }
1995  return dx;
1996 }
1997 
1998 ////////////////////////////////////////////////////////////////////////////////
1999 /// Fill vector param[4] with the bounding cylinder parameters. The order
2000 /// is the following : Rmin, Rmax, Phi1, Phi2
2001 
2003 {
2004  param[0] = TMath::Min(fRmin1, fRmin2); // Rmin
2005  param[0] *= param[0];
2006  param[1] = TMath::Max(fRmax1, fRmax2); // Rmax
2007  param[1] *= param[1];
2008  param[2] = (fPhi1<0)?(fPhi1+360.):fPhi1; // Phi1
2009  param[3] = fPhi2; // Phi2
2010  while (param[3]<param[2]) param[3]+=360.;
2011 }
2012 
2013 ////////////////////////////////////////////////////////////////////////////////
2014 /// in case shape has some negative parameters, these has to be computed
2015 /// in order to fit the mother
2016 
2018 {
2019  if (!TestShapeBit(kGeoRunTimeShape)) return 0;
2020  if (!mother->TestShapeBit(kGeoConeSeg)) {
2021  Error("GetMakeRuntimeShape", "invalid mother");
2022  return 0;
2023  }
2024  Double_t rmin1, rmax1, rmin2, rmax2, dz;
2025  rmin1 = fRmin1;
2026  rmax1 = fRmax1;
2027  rmin2 = fRmin2;
2028  rmax2 = fRmax2;
2029  dz = fDz;
2030  if (fDz<0) dz=((TGeoCone*)mother)->GetDz();
2031  if (fRmin1<0)
2032  rmin1 = ((TGeoCone*)mother)->GetRmin1();
2033  if ((fRmax1<0) || (fRmax1<fRmin1))
2034  rmax1 = ((TGeoCone*)mother)->GetRmax1();
2035  if (fRmin2<0)
2036  rmin2 = ((TGeoCone*)mother)->GetRmin2();
2037  if ((fRmax2<0) || (fRmax2<fRmin2))
2038  rmax2 = ((TGeoCone*)mother)->GetRmax2();
2039 
2040  return (new TGeoConeSeg(GetName(), dz, rmin1, rmax1, rmin2, rmax2, fPhi1, fPhi2));
2041 }
2042 
2043 ////////////////////////////////////////////////////////////////////////////////
2044 /// print shape parameters
2045 
2047 {
2048  printf("*** Shape %s: TGeoConeSeg ***\n", GetName());
2049  printf(" dz = %11.5f\n", fDz);
2050  printf(" Rmin1 = %11.5f\n", fRmin1);
2051  printf(" Rmax1 = %11.5f\n", fRmax1);
2052  printf(" Rmin2 = %11.5f\n", fRmin2);
2053  printf(" Rmax2 = %11.5f\n", fRmax2);
2054  printf(" phi1 = %11.5f\n", fPhi1);
2055  printf(" phi2 = %11.5f\n", fPhi2);
2056  printf(" Bounding box:\n");
2058 }
2059 
2060  ///////////////////////////////////////////////////////////////////////////////
2061  /// Creates a TBuffer3D describing *this* shape.
2062  /// Coordinates are in local reference frame.
2063 
2065 {
2067  Int_t nbPnts = 4*n;
2068  Int_t nbSegs = 2*nbPnts;
2069  Int_t nbPols = nbPnts-2;
2071  nbPnts, 3*nbPnts,
2072  nbSegs, 3*nbSegs,
2073  nbPols, 6*nbPols);
2074 
2075  if (buff)
2076  {
2077  SetPoints(buff->fPnts);
2078  SetSegsAndPols(*buff);
2079  }
2080 
2081  return buff;
2082 }
2083 
2084 ////////////////////////////////////////////////////////////////////////////////
2085 /// Fill TBuffer3D structure for segments and polygons.
2086 
2088 {
2089  Int_t i, j;
2091  Int_t c = GetBasicColor();
2092 
2093  memset(buffer.fSegs, 0, buffer.NbSegs()*3*sizeof(Int_t));
2094  for (i = 0; i < 4; i++) {
2095  for (j = 1; j < n; j++) {
2096  buffer.fSegs[(i*n+j-1)*3 ] = c;
2097  buffer.fSegs[(i*n+j-1)*3+1] = i*n+j-1;
2098  buffer.fSegs[(i*n+j-1)*3+2] = i*n+j;
2099  }
2100  }
2101  for (i = 4; i < 6; i++) {
2102  for (j = 0; j < n; j++) {
2103  buffer.fSegs[(i*n+j)*3 ] = c+1;
2104  buffer.fSegs[(i*n+j)*3+1] = (i-4)*n+j;
2105  buffer.fSegs[(i*n+j)*3+2] = (i-2)*n+j;
2106  }
2107  }
2108  for (i = 6; i < 8; i++) {
2109  for (j = 0; j < n; j++) {
2110  buffer.fSegs[(i*n+j)*3 ] = c;
2111  buffer.fSegs[(i*n+j)*3+1] = 2*(i-6)*n+j;
2112  buffer.fSegs[(i*n+j)*3+2] = (2*(i-6)+1)*n+j;
2113  }
2114  }
2115 
2116  Int_t indx = 0;
2117  memset(buffer.fPols, 0, buffer.NbPols()*6*sizeof(Int_t));
2118  i = 0;
2119  for (j = 0; j < n-1; j++) {
2120  buffer.fPols[indx++] = c;
2121  buffer.fPols[indx++] = 4;
2122  buffer.fPols[indx++] = (4+i)*n+j+1;
2123  buffer.fPols[indx++] = (2+i)*n+j;
2124  buffer.fPols[indx++] = (4+i)*n+j;
2125  buffer.fPols[indx++] = i*n+j;
2126  }
2127  i = 1;
2128  for (j = 0; j < n-1; j++) {
2129  buffer.fPols[indx++] = c;
2130  buffer.fPols[indx++] = 4;
2131  buffer.fPols[indx++] = i*n+j;
2132  buffer.fPols[indx++] = (4+i)*n+j;
2133  buffer.fPols[indx++] = (2+i)*n+j;
2134  buffer.fPols[indx++] = (4+i)*n+j+1;
2135  }
2136  i = 2;
2137  for (j = 0; j < n-1; j++) {
2138  buffer.fPols[indx++] = c+i;
2139  buffer.fPols[indx++] = 4;
2140  buffer.fPols[indx++] = (i-2)*2*n+j;
2141  buffer.fPols[indx++] = (4+i)*n+j;
2142  buffer.fPols[indx++] = ((i-2)*2+1)*n+j;
2143  buffer.fPols[indx++] = (4+i)*n+j+1;
2144  }
2145  i = 3;
2146  for (j = 0; j < n-1; j++) {
2147  buffer.fPols[indx++] = c+i;
2148  buffer.fPols[indx++] = 4;
2149  buffer.fPols[indx++] = (4+i)*n+j+1;
2150  buffer.fPols[indx++] = ((i-2)*2+1)*n+j;
2151  buffer.fPols[indx++] = (4+i)*n+j;
2152  buffer.fPols[indx++] = (i-2)*2*n+j;
2153  }
2154  buffer.fPols[indx++] = c+2;
2155  buffer.fPols[indx++] = 4;
2156  buffer.fPols[indx++] = 6*n;
2157  buffer.fPols[indx++] = 4*n;
2158  buffer.fPols[indx++] = 7*n;
2159  buffer.fPols[indx++] = 5*n;
2160  buffer.fPols[indx++] = c+2;
2161  buffer.fPols[indx++] = 4;
2162  buffer.fPols[indx++] = 6*n-1;
2163  buffer.fPols[indx++] = 8*n-1;
2164  buffer.fPols[indx++] = 5*n-1;
2165  buffer.fPols[indx++] = 7*n-1;
2166 }
2167 
2168 ////////////////////////////////////////////////////////////////////////////////
2169 /// computes the closest distance from given point to this shape, according
2170 /// to option. The matching point on the shape is stored in spoint.
2171 
2173 {
2174  Double_t safe = TGeoCone::Safety(point,in);
2175  if ((fPhi2-fPhi1)>=360.) return safe;
2176  Double_t safphi = TGeoShape::SafetyPhi(point, in, fPhi1, fPhi2);
2177  if (in) return TMath::Min(safe, safphi);
2178  if (safe>1.E10) return safphi;
2179  return TMath::Max(safe, safphi);
2180 }
2181 
2182 ////////////////////////////////////////////////////////////////////////////////
2183 /// Static method to compute the closest distance from given point to this shape.
2184 
2186  Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2, Int_t skipz)
2187 {
2188  Double_t safe = TGeoCone::SafetyS(point,in,dz,rmin1,rmax1,rmin2,rmax2,skipz);
2189  if ((phi2-phi1)>=360.) return safe;
2190  Double_t safphi = TGeoShape::SafetyPhi(point,in,phi1,phi2);
2191  if (in) return TMath::Min(safe, safphi);
2192  if (safe>1.E10) return safphi;
2193  return TMath::Max(safe, safphi);
2194 }
2195 
2196 ////////////////////////////////////////////////////////////////////////////////
2197 /// Save a primitive as a C++ statement(s) on output stream "out".
2198 
2199 void TGeoConeSeg::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
2200 {
2201  if (TObject::TestBit(kGeoSavePrimitive)) return;
2202  out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
2203  out << " dz = " << fDz << ";" << std::endl;
2204  out << " rmin1 = " << fRmin1 << ";" << std::endl;
2205  out << " rmax1 = " << fRmax1 << ";" << std::endl;
2206  out << " rmin2 = " << fRmin2 << ";" << std::endl;
2207  out << " rmax2 = " << fRmax2 << ";" << std::endl;
2208  out << " phi1 = " << fPhi1 << ";" << std::endl;
2209  out << " phi2 = " << fPhi2 << ";" << std::endl;
2210  out << " TGeoShape *" << GetPointerName() << " = new TGeoConeSeg(\"" << GetName() << "\", dz,rmin1,rmax1,rmin2,rmax2,phi1,phi2);" << std::endl;
2212 }
2213 
2214 ////////////////////////////////////////////////////////////////////////////////
2215 /// Set dimensions of the cone segment.
2216 
2218  Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2)
2219 {
2220  fDz = dz;
2221  fRmin1 = rmin1;
2222  fRmax1 = rmax1;
2223  fRmin2 = rmin2;
2224  fRmax2 = rmax2;
2225  fPhi1 = phi1;
2226  while (fPhi1<0) fPhi1+=360.;
2227  fPhi2 = phi2;
2228  while (fPhi2<=fPhi1) fPhi2+=360.;
2229  if (TGeoShape::IsSameWithinTolerance(fPhi1,fPhi2)) Fatal("SetConsDimensions", "In shape %s invalid phi1=%g, phi2=%g\n", GetName(), fPhi1, fPhi2);
2230  InitTrigonometry();
2231 }
2232 
2233 ////////////////////////////////////////////////////////////////////////////////
2234 /// Set dimensions of the cone segment from an array.
2235 
2237 {
2238  Double_t dz = param[0];
2239  Double_t rmin1 = param[1];
2240  Double_t rmax1 = param[2];
2241  Double_t rmin2 = param[3];
2242  Double_t rmax2 = param[4];
2243  Double_t phi1 = param[5];
2244  Double_t phi2 = param[6];
2245  SetConsDimensions(dz, rmin1, rmax1,rmin2, rmax2, phi1, phi2);
2246 }
2247 
2248 ////////////////////////////////////////////////////////////////////////////////
2249 /// Create cone segment mesh points.
2250 
2252 {
2253  Int_t j, n;
2254  Float_t dphi,phi,phi1, phi2,dz;
2255 
2256  n = gGeoManager->GetNsegments()+1;
2257  dz = fDz;
2258  phi1 = fPhi1;
2259  phi2 = fPhi2;
2260 
2261  dphi = (phi2-phi1)/(n-1);
2262 
2263  Int_t indx = 0;
2264 
2265  if (points) {
2266  for (j = 0; j < n; j++) {
2267  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2268  points[indx++] = fRmin1 * TMath::Cos(phi);
2269  points[indx++] = fRmin1 * TMath::Sin(phi);
2270  points[indx++] = -dz;
2271  }
2272  for (j = 0; j < n; j++) {
2273  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2274  points[indx++] = fRmax1 * TMath::Cos(phi);
2275  points[indx++] = fRmax1 * TMath::Sin(phi);
2276  points[indx++] = -dz;
2277  }
2278  for (j = 0; j < n; j++) {
2279  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2280  points[indx++] = fRmin2 * TMath::Cos(phi);
2281  points[indx++] = fRmin2 * TMath::Sin(phi);
2282  points[indx++] = dz;
2283  }
2284  for (j = 0; j < n; j++) {
2285  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2286  points[indx++] = fRmax2 * TMath::Cos(phi);
2287  points[indx++] = fRmax2 * TMath::Sin(phi);
2288  points[indx++] = dz;
2289  }
2290  }
2291 }
2292 
2293 ////////////////////////////////////////////////////////////////////////////////
2294 /// Create cone segment mesh points.
2295 
2297 {
2298  Int_t j, n;
2299  Float_t dphi,phi,phi1, phi2,dz;
2300 
2301  n = gGeoManager->GetNsegments()+1;
2302  dz = fDz;
2303  phi1 = fPhi1;
2304  phi2 = fPhi2;
2305 
2306  dphi = (phi2-phi1)/(n-1);
2307 
2308  Int_t indx = 0;
2309 
2310  if (points) {
2311  for (j = 0; j < n; j++) {
2312  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2313  points[indx++] = fRmin1 * TMath::Cos(phi);
2314  points[indx++] = fRmin1 * TMath::Sin(phi);
2315  points[indx++] = -dz;
2316  }
2317  for (j = 0; j < n; j++) {
2318  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2319  points[indx++] = fRmax1 * TMath::Cos(phi);
2320  points[indx++] = fRmax1 * TMath::Sin(phi);
2321  points[indx++] = -dz;
2322  }
2323  for (j = 0; j < n; j++) {
2324  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2325  points[indx++] = fRmin2 * TMath::Cos(phi);
2326  points[indx++] = fRmin2 * TMath::Sin(phi);
2327  points[indx++] = dz;
2328  }
2329  for (j = 0; j < n; j++) {
2330  phi = (fPhi1+j*dphi)*TMath::DegToRad();
2331  points[indx++] = fRmax2 * TMath::Cos(phi);
2332  points[indx++] = fRmax2 * TMath::Sin(phi);
2333  points[indx++] = dz;
2334  }
2335  }
2336 }
2337 
2338 ////////////////////////////////////////////////////////////////////////////////
2339 /// Returns numbers of vertices, segments and polygons composing the shape mesh.
2340 
2341 void TGeoConeSeg::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
2342 {
2344  nvert = n*4;
2345  nsegs = n*8;
2346  npols = n*4-2;
2347 }
2348 
2349 ////////////////////////////////////////////////////////////////////////////////
2350 /// Return number of vertices of the mesh representation
2351 
2353 {
2355  Int_t numPoints = n*4;
2356  return numPoints;
2357 }
2358 
2359 ////////////////////////////////////////////////////////////////////////////////
2360 /// Fill size of this 3-D object
2361 
2363 {
2364 }
2365 
2366 ////////////////////////////////////////////////////////////////////////////////
2367 /// Fills a static 3D buffer and returns a reference.
2368 
2369 const TBuffer3D & TGeoConeSeg::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
2370 {
2371  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
2372 
2373  TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
2374 
2375  if (reqSections & TBuffer3D::kRawSizes) {
2377  Int_t nbPnts = 4*n;
2378  Int_t nbSegs = 2*nbPnts;
2379  Int_t nbPols = nbPnts-2;
2380  if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
2382  }
2383  }
2384  if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
2385  SetPoints(buffer.fPnts);
2386  if (!buffer.fLocalFrame) {
2387  TransformPoints(buffer.fPnts, buffer.NbPnts());
2388  }
2389 
2390  SetSegsAndPols(buffer);
2392  }
2393 
2394  return buffer;
2395 }
2396 
2397 ////////////////////////////////////////////////////////////////////////////////
2398 /// Fills array with n random points located on the line segments of the shape mesh.
2399 /// The output array must be provided with a length of minimum 3*npoints. Returns
2400 /// true if operation is implemented.
2401 
2403 {
2404  if (npoints > (npoints/2)*2) {
2405  Error("GetPointsOnSegments","Npoints must be even number");
2406  return kFALSE;
2407  }
2408  Int_t nc = (Int_t)TMath::Sqrt(0.5*npoints);
2409  Double_t dphi = (fPhi2-fPhi1)*TMath::DegToRad()/(nc-1);
2410  Double_t phi = 0;
2411  Double_t phi1 = fPhi1 * TMath::DegToRad();
2412  Int_t ntop = npoints/2 - nc*(nc-1);
2413  Double_t dz = 2*fDz/(nc-1);
2414  Double_t z = 0;
2415  Double_t rmin = 0.;
2416  Double_t rmax = 0.;
2417  Int_t icrt = 0;
2418  Int_t nphi = nc;
2419  // loop z sections
2420  for (Int_t i=0; i<nc; i++) {
2421  if (i == (nc-1)) {
2422  nphi = ntop;
2423  dphi = (fPhi2-fPhi1)*TMath::DegToRad()/(nphi-1);
2424  }
2425  z = -fDz + i*dz;
2426  rmin = 0.5*(fRmin1+fRmin2) + 0.5*(fRmin2-fRmin1)*z/fDz;
2427  rmax = 0.5*(fRmax1+fRmax2) + 0.5*(fRmax2-fRmax1)*z/fDz;
2428  // loop points on circle sections
2429  for (Int_t j=0; j<nphi; j++) {
2430  phi = phi1 + j*dphi;
2431  array[icrt++] = rmin * TMath::Cos(phi);
2432  array[icrt++] = rmin * TMath::Sin(phi);
2433  array[icrt++] = z;
2434  array[icrt++] = rmax * TMath::Cos(phi);
2435  array[icrt++] = rmax * TMath::Sin(phi);
2436  array[icrt++] = z;
2437  }
2438  }
2439  return kTRUE;
2440 }
2441 
2442 ////////////////////////////////////////////////////////////////////////////////
2443 /// Check the inside status for each of the points in the array.
2444 /// Input: Array of point coordinates + vector size
2445 /// Output: Array of Booleans for the inside of each point
2446 
2447 void TGeoConeSeg::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
2448 {
2449  for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
2450 }
2451 
2452 ////////////////////////////////////////////////////////////////////////////////
2453 /// Compute the normal for an array o points so that norm.dot.dir is positive
2454 /// Input: Arrays of point coordinates and directions + vector size
2455 /// Output: Array of normal directions
2456 
2457 void TGeoConeSeg::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
2458 {
2459  for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
2460 }
2461 
2462 ////////////////////////////////////////////////////////////////////////////////
2463 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
2464 
2465 void TGeoConeSeg::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2466 {
2467  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
2468 }
2469 
2470 ////////////////////////////////////////////////////////////////////////////////
2471 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
2472 
2473 void TGeoConeSeg::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2474 {
2475  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
2476 }
2477 
2478 ////////////////////////////////////////////////////////////////////////////////
2479 /// Compute safe distance from each of the points in the input array.
2480 /// Input: Array of point coordinates, array of statuses for these points, size of the arrays
2481 /// Output: Safety values
2482 
2483 void TGeoConeSeg::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
2484 {
2485  for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
2486 }
c
#define c(i)
Definition: RSha256.hxx:119
TGeoConeSeg::Contains_v
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: TGeoCone.cxx:2447
n
const Int_t n
Definition: legend1.C:16
TGeoCone::SetSegsAndPols
virtual void SetSegsAndPols(TBuffer3D &buffer) const
Fill TBuffer3D structure for segments and polygons.
Definition: TGeoCone.cxx:787
TGeoCone::Capacity
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoCone.cxx:145
TGeoCone::GetBoundingCylinder
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoCone.cxx:657
TGeoConeSeg::GetBuffer3D
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
Definition: TGeoCone.cxx:2369
ymax
float ymax
Definition: THbookFile.cxx:95
TGeoConeSeg::fSm
Double_t fSm
Definition: TGeoCone.h:109
TGeoCone::GetAxisRange
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
Definition: TGeoCone.cxx:634
TGeoVolume::SetFinder
void SetFinder(TGeoPatternFinder *finder)
Definition: TGeoVolume.h:231
TGeoConeSeg::GetBoundingCylinder
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoCone.cxx:2002
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
TGeoConeSeg::~TGeoConeSeg
virtual ~TGeoConeSeg()
destructor
Definition: TGeoCone.cxx:1256
TGeoConeSeg::Capacity
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoCone.cxx:1289
TBuffer3D::SectionsValid
Bool_t SectionsValid(UInt_t mask) const
Definition: TBuffer3D.h:73
TMath::ATan2
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:679
TMath::RadToDeg
constexpr Double_t RadToDeg()
Conversion from radian to degree:
Definition: TMath.h:79
TGeoVolume::GetNdaughters
Int_t GetNdaughters() const
Definition: TGeoVolume.h:349
TBuffer3D::SetSectionsValid
void SetSectionsValid(UInt_t mask)
Definition: TBuffer3D.h:71
TMath::Max
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
TGeoConeSeg::Sizeof3D
virtual void Sizeof3D() const
Fill size of this 3-D object.
Definition: TGeoCone.cxx:2362
TGeoConeSeg::Divide
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this cone segment shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition: TGeoCone.cxx:1918
TMath::Cos
Double_t Cos(Double_t)
Definition: TMath.h:643
gGeoManager
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
TString::Data
const char * Data() const
Definition: TString.h:369
TGeoConeSeg::GetNmeshVertices
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
Definition: TGeoCone.cxx:2352
TGeoCone::GetMeshNumbers
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: TGeoCone.cxx:1090
TGeoConeSeg::DistFromOutsideS
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Double_t cm, Double_t sm, Double_t cdfi)
compute distance from outside point to surface of arbitrary tube
Definition: TGeoCone.cxx:1610
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TGeoShape::SafetyPhi
static Double_t SafetyPhi(const Double_t *point, Bool_t in, Double_t phi1, Double_t phi2)
Static method to compute safety w.r.t a phi corner defined by cosines/sines of the angles phi1,...
Definition: TGeoShape.cxx:464
TGeoCone::Safety
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: TGeoCone.cxx:867
TGeoPatternCylPhi
Definition: TGeoPatternFinder.h:395
r
ROOT::R::TRInterface & r
Definition: Object.C:4
xmax
float xmax
Definition: THbookFile.cxx:95
TGeoConeSeg::ComputeBBox
virtual void ComputeBBox()
compute bounding box of the tube segment
Definition: TGeoCone.cxx:1308
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:890
TGeoCone::DistFromOutside_v
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: TGeoCone.cxx:1179
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
TGeoConeSeg::fS1
Double_t fS1
Definition: TGeoCone.h:105
TGeoConeSeg::Safety
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: TGeoCone.cxx:2172
TMath::DegToRad
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:87
TGeoBBox::fOrigin
Double_t fOrigin[3]
Definition: TGeoBBox.h:30
Float_t
float Float_t
Definition: RtypesCore.h:57
TGeoConeSeg::fC2
Double_t fC2
Definition: TGeoCone.h:108
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:168
TGeoShape::TransformPoints
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
Definition: TGeoShape.cxx:552
TObject::Fatal
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition: TObject.cxx:918
Int_t
int Int_t
Definition: RtypesCore.h:45
TGeoVolume.h
box
void box(Int_t pat, Double_t x1, Double_t y1, Double_t x2, Double_t y2)
Definition: fillpatterns.C:1
TGeoConeSeg::InspectShape
virtual void InspectShape() const
print shape parameters
Definition: TGeoCone.cxx:2046
TGeoManager::GetNsegments
Int_t GetNsegments() const
Get number of segments approximating circles.
Definition: TGeoManager.cxx:3349
TGeoConeSeg::InitTrigonometry
void InitTrigonometry()
Init frequently used trigonometric values.
Definition: TGeoCone.cxx:1271
TBuffer3D::NbPnts
UInt_t NbPnts() const
Definition: TBuffer3D.h:86
TGeoCone::MakeBuffer3D
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
Definition: TGeoCone.cxx:765
TGeoConeSeg::Contains
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this sphere
Definition: TGeoCone.cxx:1452
TGeoNodeOffset
Definition: TGeoNode.h:183
TGeoConeSeg::SetSegsAndPols
virtual void SetSegsAndPols(TBuffer3D &buffer) const
Fill TBuffer3D structure for segments and polygons.
Definition: TGeoCone.cxx:2087
TGeoShape::SetShapeBit
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
Definition: TGeoShape.cxx:524
TGeoCone
Definition: TGeoCone.h:17
TBuffer3D::SetRawSizes
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
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TGeoShape::GetBasicColor
Int_t GetBasicColor() const
Get the basic color (0-7).
Definition: TGeoShape.cxx:673
TGeoCone::InspectShape
virtual void InspectShape() const
print shape parameters
Definition: TGeoCone.cxx:749
TGeoConeSeg::DistFromInside_v
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: TGeoCone.cxx:2465
TObjArray::At
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
TBuffer3D::fSegs
Int_t * fSegs
Definition: TBuffer3D.h:119
TString
Definition: TString.h:136
TMath::LocMax
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition: TMath.h:1000
TGeoCone::Contains_v
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: TGeoCone.cxx:1153
TGeoCone.h
TGeoCone::ComputeNormal
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoCone.cxx:180
TGeant4Unit::cm
static constexpr double cm
Definition: TGeant4SystemOfUnits.h:118
TGeoConeSeg::GetMeshNumbers
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: TGeoCone.cxx:2341
b
#define b(i)
Definition: RSha256.hxx:118
TGeoConeSeg::fPhi1
Double_t fPhi1
Definition: TGeoCone.h:102
TGeoPatternFinder
Definition: TGeoPatternFinder.h:31
bool
TGeoShape::TestShapeBit
Bool_t TestShapeBit(UInt_t f) const
Definition: TGeoShape.h:163
TGeoConeSeg::Safety_v
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: TGeoCone.cxx:2483
TGeoVolume::GetMedium
TGeoMedium * GetMedium() const
Definition: TGeoVolume.h:173
TGeoCone::DistFromInside
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 cone Boundary safe algorithm.
Definition: TGeoCone.cxx:346
TGeoShape::NormalPhi
static void NormalPhi(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Static method to compute normal to phi planes.
Definition: TGeoShape.cxx:437
id
XFontStruct * id
Definition: TGX11.cxx:109
TGeoConeSeg::SafetyS
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2, Int_t skipz=0)
Static method to compute the closest distance from given point to this shape.
Definition: TGeoCone.cxx:2185
TGeoVolumeMulti::AddVolume
void AddVolume(TGeoVolume *vol)
Add a volume with valid shape to the list of volumes.
Definition: TGeoVolume.cxx:2419
TGeoCone::fRmax2
Double_t fRmax2
Definition: TGeoCone.h:31
TGeoCone::DistFromOutsideS
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
Compute distance from outside point to surface of the tube Boundary safe algorithm.
Definition: TGeoCone.cxx:361
TGeoShape::SafetySeg
static Double_t SafetySeg(Double_t r, Double_t z, Double_t r1, Double_t z1, Double_t r2, Double_t z2, Bool_t outer)
Compute distance from point of coordinates (r,z) to segment (r1,z1):(r2,z2)
Definition: TGeoShape.cxx:494
TGeoConeSeg::fS2
Double_t fS2
Definition: TGeoCone.h:107
TGeoCone::fDz
Double_t fDz
Definition: TGeoCone.h:27
TGeoConeSeg::fC1
Double_t fC1
Definition: TGeoCone.h:106
TGeoConeSeg::TGeoConeSeg
TGeoConeSeg()
Default constructor.
Definition: TGeoCone.cxx:1199
TGeoShape::IsCloseToPhi
static Bool_t IsCloseToPhi(Double_t epsil, const Double_t *point, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
True if point is closer than epsil to one of the phi planes defined by c1,s1 or c2,...
Definition: TGeoShape.cxx:269
TMath::LocMin
Long64_t LocMin(Long64_t n, const T *a)
Return index of array with the minimum element.
Definition: TMath.h:972
TBuffer3D
Definition: TBuffer3D.h:17
TMath::Pi
constexpr Double_t Pi()
Definition: TMath.h:43
TGeoConeSeg::AfterStreamer
virtual void AfterStreamer()
Function called after streaming an object of this class.
Definition: TGeoCone.cxx:1263
Option_t
const typedef char Option_t
Definition: RtypesCore.h:66
TGeoShape
Definition: TGeoShape.h:25
TGeoConeSeg::GetMakeRuntimeShape
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: TGeoCone.cxx:2017
TGeoCone::DistancetoPrimitive
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each corner
Definition: TGeoCone.cxx:545
TGeoCone::SetPoints
virtual void SetPoints(Double_t *points) const
Create cone mesh points.
Definition: TGeoCone.cxx:1002
xmin
float xmin
Definition: THbookFile.cxx:95
TBuffer3DTypes::kGeneric
@ kGeneric
Definition: TBuffer3DTypes.h:36
TGeoBBox::InspectShape
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoBBox.cxx:792
TGeoCone::GetBuffer3D
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
Definition: TGeoCone.cxx:1118
TObject::SetBit
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
TGeoBBox::fDY
Double_t fDY
Definition: TGeoBBox.h:28
TGeoVolume::GetNodes
TObjArray * GetNodes()
Definition: TGeoVolume.h:167
a
auto * a
Definition: textangle.C:12
TGeoCone::SafetyS
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Int_t skipz=0)
computes the closest distance from given point to this shape, according to option.
Definition: TGeoCone.cxx:884
TGeoCone::Sizeof3D
virtual void Sizeof3D() const
Fill size of this 3-D object.
Definition: TGeoCone.cxx:1111
TBuffer3D.h
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TMath::Sign
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:165
s1
#define s1(x)
Definition: RSha256.hxx:109
TGeoShape::kGeoSavePrimitive
@ kGeoSavePrimitive
Definition: TGeoShape.h:65
TGeant4Unit::sr
static constexpr double sr
Definition: TGeant4SystemOfUnits.h:150
TGeoConeSeg::SetDimensions
virtual void SetDimensions(Double_t *param)
Set dimensions of the cone segment from an array.
Definition: TGeoCone.cxx:2236
TGeoCone::GetAxisName
virtual const char * GetAxisName(Int_t iaxis) const
Returns name of axis IAXIS.
Definition: TGeoCone.cxx:617
TBuffer3DTypes.h
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:639
TGeoShape::kGeoRunTimeShape
@ kGeoRunTimeShape
Definition: TGeoShape.h:41
TBuffer3D::kRaw
@ kRaw
Definition: TBuffer3D.h:60
TGeoPatternFinder::SetDivIndex
void SetDivIndex(Int_t index)
Definition: TGeoPatternFinder.h:99
TGeoCone::fRmin2
Double_t fRmin2
Definition: TGeoCone.h:30
TGeoShape::GetName
virtual const char * GetName() const
Get the shape name.
Definition: TGeoShape.cxx:248
TGeoBBox
Definition: TGeoBBox.h:17
TGeoConeSeg::SetPoints
virtual void SetPoints(Double_t *points) const
Create cone segment mesh points.
Definition: TGeoCone.cxx:2251
TGeoCone::Contains
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this cone
Definition: TGeoCone.cxx:260
TGeoBBox::FillBuffer3D
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
TObject::Warning
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:876
TMath::Min
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
TGeoCone::ComputeNormal_v
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: TGeoCone.cxx:1163
TGeoCone::GetPointsOnSegments
virtual Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const
Fills array with n random points located on the line segments of the shape mesh.
Definition: TGeoCone.cxx:702
TMath::TwoPi
constexpr Double_t TwoPi()
Definition: TMath.h:50
ymin
float ymin
Definition: THbookFile.cxx:95
TGeoConeSeg::MakeBuffer3D
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
Definition: TGeoCone.cxx:2064
TGeoConeSeg::DistancetoPrimitive
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each corner
Definition: TGeoCone.cxx:1903
TBuffer3D::NbPols
UInt_t NbPols() const
Definition: TBuffer3D.h:88
TGeoManager.h
TGeoBBox::fDZ
Double_t fDZ
Definition: TGeoBBox.h:29
TGeoConeSeg::DistToCons
static Double_t DistToCons(const Double_t *point, const Double_t *dir, Double_t r1, Double_t z1, Double_t r2, Double_t z2, Double_t phi1, Double_t phi2)
Static method to compute distance to a conical surface with :
Definition: TGeoCone.cxx:1472
TGeoCone::DistFromOutside
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
Definition: TGeoCone.cxx:491
TGeoShape::GetPointerName
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:699
TGeoCone::DistFromInsideS
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
Compute distance from inside point to surface of the cone (static) Boundary safe algorithm.
Definition: TGeoCone.cxx:274
TVirtualGeoPainter.h
TGeoCone::GetMakeRuntimeShape
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: TGeoCone.cxx:671
Double_t
double Double_t
Definition: RtypesCore.h:59
TGeoConeSeg::DistFromInsideS
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Double_t cm, Double_t sm, Double_t cdfi)
compute distance from inside point to surface of the tube segment
Definition: TGeoCone.cxx:1544
TGeoMatrix
Definition: TGeoMatrix.h:40
TBuffer3D::fLocalFrame
Bool_t fLocalFrame
Definition: TBuffer3D.h:96
TGeoConeSeg::ComputeNormal
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoCone.cxx:1354
TBuffer3D::NbSegs
UInt_t NbSegs() const
Definition: TBuffer3D.h:87
TGeoConeSeg::SetConsDimensions
void SetConsDimensions(Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t phi1, Double_t phi2)
Set dimensions of the cone segment.
Definition: TGeoCone.cxx:2217
TGeoManager::MakeVolumeMulti
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Definition: TGeoManager.cxx:3308
TGeoCone::SetDimensions
virtual void SetDimensions(Double_t *param)
Set cone dimensions from an array.
Definition: TGeoCone.cxx:989
TGeoCone::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoCone.cxx:919
TGeoVolumeMulti
Definition: TGeoVolume.h:253
TGeoCone::ComputeBBox
virtual void ComputeBBox()
compute bounding box of the sphere
Definition: TGeoCone.cxx:170
TGeoCone::Safety_v
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: TGeoCone.cxx:1189
TGeoPatternZ
Definition: TGeoPatternFinder.h:184
TGeoConeSeg::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoCone.cxx:2199
points
point * points
Definition: X3DBuffer.c:22
TGeoShape::IsSameWithinTolerance
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
Definition: TGeoShape.cxx:326
TGeoCone::DistFromInside_v
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: TGeoCone.cxx:1171
TGeoCone::fRmin1
Double_t fRmin1
Definition: TGeoCone.h:28
TGeoCone::~TGeoCone
virtual ~TGeoCone()
destructor
Definition: TGeoCone.cxx:163
TGeoCone::fRmax1
Double_t fRmax1
Definition: TGeoCone.h:29
TBuffer3D::kRawSizes
@ kRawSizes
Definition: TBuffer3D.h:59
name
char name[80]
Definition: TGX11.cxx:110
TGeoConeSeg::ComputeNormalS
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Compute normal to closest surface from POINT.
Definition: TGeoCone.cxx:1405
TGeoShape::ShapeDistancetoPrimitive
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
Definition: TGeoShape.cxx:259
TGeoConeSeg::ComputeNormal_v
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: TGeoCone.cxx:2457
TGeoShape::DistToPhiMin
static Double_t DistToPhiMin(const Double_t *point, const Double_t *dir, Double_t s1, Double_t c1, Double_t s2, Double_t c2, Double_t sm, Double_t cm, Bool_t in=kTRUE)
compute distance from point (inside phi) to both phi planes. Return minimum.
Definition: TGeoShape.cxx:405
c2
return c2
Definition: legend2.C:14
TGeoShape::kGeoBad
@ kGeoBad
Definition: TGeoShape.h:34
TBuffer3D::fPnts
Double_t * fPnts
Definition: TBuffer3D.h:118
TGeoShape::kGeoConeSeg
@ kGeoConeSeg
Definition: TGeoShape.h:50
TGeoConeSeg
Definition: TGeoCone.h:98
TGeoConeSeg::DistFromOutside
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
Definition: TGeoCone.cxx:1885
TGeoCone::SetConeDimensions
void SetConeDimensions(Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
Set cone dimensions.
Definition: TGeoCone.cxx:935
TGeoConeSeg::fCdfi
Double_t fCdfi
Definition: TGeoCone.h:111
TGeoCone::DistToCone
static void DistToCone(const Double_t *point, const Double_t *dir, Double_t dz, Double_t r1, Double_t r2, Double_t &b, Double_t &delta)
Static method to compute distance to a conical surface with :
Definition: TGeoCone.cxx:511
TGeoConeSeg::fPhi2
Double_t fPhi2
Definition: TGeoCone.h:103
TGeoConeSeg::DistFromInside
virtual Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=0) const
compute distance from inside point to surface of the tube segment
Definition: TGeoCone.cxx:1594
TGeoShape::Tolerance
static Double_t Tolerance()
Definition: TGeoShape.h:91
TGeoCone::ComputeNormalS
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t dz, Double_t rmin1, Double_t rmax1, Double_t rmin2, Double_t rmax2)
Compute normal to closest surface from POINT.
Definition: TGeoCone.cxx:223
TGeoConeSeg::DistFromOutside_v
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: TGeoCone.cxx:2473
TGeoConeSeg::GetAxisRange
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
Definition: TGeoCone.cxx:1978
TGeoConeSeg::fCm
Double_t fCm
Definition: TGeoCone.h:110
TGeoShape::Big
static Double_t Big()
Definition: TGeoShape.h:88
TGeoVolume::AddNodeOffset
void AddNodeOffset(TGeoVolume *vol, Int_t copy_no, Double_t offset=0, Option_t *option="")
Add a division node to the list of nodes.
Definition: TGeoVolume.cxx:971
TBuffer3D::fPols
Int_t * fPols
Definition: TBuffer3D.h:120
TObject::ClassName
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
TGeoBBox::DistFromOutside
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
TGeoShape::kGeoCone
@ kGeoCone
Definition: TGeoShape.h:49
TGeoVolume
Definition: TGeoVolume.h:44
TGeoBBox::fDX
Double_t fDX
Definition: TGeoBBox.h:27
TGeoCone::GetNmeshVertices
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
Definition: TGeoCone.cxx:1101
TMath.h
TGeoConeSeg::GetPointsOnSegments
virtual Bool_t GetPointsOnSegments(Int_t npoints, Double_t *array) const
Fills array with n random points located on the line segments of the shape mesh.
Definition: TGeoCone.cxx:2402
int
c1
return c1
Definition: legend1.C:41
TGeoCone::Divide
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this cone shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition: TGeoCone.cxx:560
TGeoCone::TGeoCone
TGeoCone()
Default constructor.
Definition: TGeoCone.cxx:84