Logo ROOT  
Reference Guide
TGeoPcon.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 24/10/01
3 // TGeoPcon::Contains() 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 TGeoPcon
14 \ingroup Geometry_classes
15 A polycone.
16 
17 It has at least 9 parameters:
18  - the lower phi limit;
19  - the range in phi;
20  - the number of z planes (at least two) where the inner/outer
21  radii are changing;
22  - z coordinate, inner and outer radius for each z plane
23 
24 Begin_Macro(source)
25 {
26  TCanvas *c = new TCanvas("c", "c",0,0,600,600);
27  new TGeoManager("pcon", "poza10");
28  TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
29  TGeoMedium *med = new TGeoMedium("MED",1,mat);
30  TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
31  gGeoManager->SetTopVolume(top);
32  TGeoVolume *vol = gGeoManager->MakePcon("PCON",med, -30.0,300,4);
33  TGeoPcon *pcon = (TGeoPcon*)(vol->GetShape());
34  pcon->DefineSection(0,0,15,20);
35  pcon->DefineSection(1,20,15,20);
36  pcon->DefineSection(2,20,15,25);
37  pcon->DefineSection(3,50,15,20);
38  vol->SetLineWidth(2);
39  top->AddNode(vol,1);
40  gGeoManager->CloseGeometry();
41  gGeoManager->SetNsegments(30);
42  top->Draw();
43  TView *view = gPad->GetView();
44  view->ShowAxis();
45 }
46 End_Macro
47 */
48 
49 #include "TGeoPcon.h"
50 
51 #include <iostream>
52 
53 #include "TBuffer.h"
54 #include "TGeoManager.h"
55 #include "TGeoVolume.h"
56 #include "TVirtualGeoPainter.h"
57 #include "TGeoTube.h"
58 #include "TGeoCone.h"
59 #include "TBuffer3D.h"
60 #include "TBuffer3DTypes.h"
61 #include "TMath.h"
62 
64 
65 ////////////////////////////////////////////////////////////////////////////////
66 /// dummy ctor
67 
69  :TGeoBBox(),
70  fNz(0),
71  fPhi1(0.),
72  fDphi(0.),
73  fRmin(nullptr),
74  fRmax(nullptr),
75  fZ(nullptr),
76  fFullPhi(kFALSE),
77  fC1(0.),
78  fS1(0.),
79  fC2(0.),
80  fS2(0.),
81  fCm(0.),
82  fSm(0.),
83  fCdphi(0.)
84 {
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 /// Default constructor
90 
92  :TGeoBBox(0, 0, 0),
93  fNz(nz),
94  fPhi1(phi),
95  fDphi(dphi),
96  fRmin(nullptr),
97  fRmax(nullptr),
98  fZ(nullptr),
99  fFullPhi(kFALSE),
100  fC1(0.),
101  fS1(0.),
102  fC2(0.),
103  fS2(0.),
104  fCm(0.),
105  fSm(0.),
106  fCdphi(0.)
107 {
109  while (fPhi1<0) fPhi1+=360.;
110  fRmin = new Double_t [nz];
111  fRmax = new Double_t [nz];
112  fZ = new Double_t [nz];
113  memset(fRmin, 0, nz*sizeof(Double_t));
114  memset(fRmax, 0, nz*sizeof(Double_t));
115  memset(fZ, 0, nz*sizeof(Double_t));
117  Double_t phi1 = fPhi1;
118  Double_t phi2 = phi1+fDphi;
119  Double_t phim = 0.5*(phi1+phi2);
120  fC1 = TMath::Cos(phi1*TMath::DegToRad());
121  fS1 = TMath::Sin(phi1*TMath::DegToRad());
122  fC2 = TMath::Cos(phi2*TMath::DegToRad());
123  fS2 = TMath::Sin(phi2*TMath::DegToRad());
124  fCm = TMath::Cos(phim*TMath::DegToRad());
125  fSm = TMath::Sin(phim*TMath::DegToRad());
127 }
128 
129 ////////////////////////////////////////////////////////////////////////////////
130 /// Default constructor
131 
132 TGeoPcon::TGeoPcon(const char *name, Double_t phi, Double_t dphi, Int_t nz)
133  :TGeoBBox(name, 0, 0, 0),
134  fNz(nz),
135  fPhi1(phi),
136  fDphi(dphi),
137  fRmin(nullptr),
138  fRmax(nullptr),
139  fZ(nullptr),
140  fFullPhi(kFALSE),
141  fC1(0.),
142  fS1(0.),
143  fC2(0.),
144  fS2(0.),
145  fCm(0.),
146  fSm(0.),
147  fCdphi(0.)
148 {
150  while (fPhi1<0) fPhi1+=360.;
151  fRmin = new Double_t [nz];
152  fRmax = new Double_t [nz];
153  fZ = new Double_t [nz];
154  memset(fRmin, 0, nz*sizeof(Double_t));
155  memset(fRmax, 0, nz*sizeof(Double_t));
156  memset(fZ, 0, nz*sizeof(Double_t));
158  Double_t phi1 = fPhi1;
159  Double_t phi2 = phi1+fDphi;
160  Double_t phim = 0.5*(phi1+phi2);
161  fC1 = TMath::Cos(phi1*TMath::DegToRad());
162  fS1 = TMath::Sin(phi1*TMath::DegToRad());
163  fC2 = TMath::Cos(phi2*TMath::DegToRad());
164  fS2 = TMath::Sin(phi2*TMath::DegToRad());
165  fCm = TMath::Cos(phim*TMath::DegToRad());
166  fSm = TMath::Sin(phim*TMath::DegToRad());
168 }
169 
170 ////////////////////////////////////////////////////////////////////////////////
171 /// Default constructor in GEANT3 style
172 /// - param[0] = phi1
173 /// - param[1] = dphi
174 /// - param[2] = nz
175 /// - param[3] = z1
176 /// - param[4] = Rmin1
177 /// - param[5] = Rmax1
178 /// ...
179 
181  :TGeoBBox(0, 0, 0),
182  fNz(0),
183  fPhi1(0.),
184  fDphi(0.),
185  fRmin(nullptr),
186  fRmax(nullptr),
187  fZ(nullptr),
188  fFullPhi(kFALSE),
189  fC1(0.),
190  fS1(0.),
191  fC2(0.),
192  fS2(0.),
193  fCm(0.),
194  fSm(0.),
195  fCdphi(0.)
196 {
198  SetDimensions(param);
199  ComputeBBox();
200 }
201 
202 ////////////////////////////////////////////////////////////////////////////////
203 /// destructor
204 
206 {
207  if (fRmin) { delete[] fRmin; fRmin = nullptr; }
208  if (fRmax) { delete[] fRmax; fRmax = nullptr; }
209  if (fZ) { delete[] fZ; fZ = nullptr; }
210 }
211 
212 ////////////////////////////////////////////////////////////////////////////////
213 /// Computes capacity of the shape in [length^3]
214 
216 {
217  Int_t ipl;
218  Double_t rmin1, rmax1, rmin2, rmax2, phi1, phi2, dz;
219  Double_t capacity = 0.;
220  phi1 = fPhi1;
221  phi2 = fPhi1 + fDphi;
222  for (ipl=0; ipl<fNz-1; ipl++) {
223  dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
224  if (dz < TGeoShape::Tolerance()) continue;
225  rmin1 = fRmin[ipl];
226  rmax1 = fRmax[ipl];
227  rmin2 = fRmin[ipl+1];
228  rmax2 = fRmax[ipl+1];
229  capacity += TGeoConeSeg::Capacity(dz,rmin1,rmax1,rmin2,rmax2,phi1,phi2);
230  }
231  return capacity;
232 }
233 
234 ////////////////////////////////////////////////////////////////////////////////
235 /// compute bounding box of the pcon
236 /// Check if the sections are in increasing Z order
237 
239 {
240  for (Int_t isec=0; isec<fNz-1; isec++) {
241  if (TMath::Abs(fZ[isec]-fZ[isec+1]) < TGeoShape::Tolerance()) {
242  fZ[isec+1]=fZ[isec];
243  if (IsSameWithinTolerance(fRmin[isec], fRmin[isec+1]) &&
244  IsSameWithinTolerance(fRmax[isec], fRmax[isec+1])) {
245  InspectShape();
246  Error("ComputeBBox", "Duplicated section %d/%d for shape %s", isec, isec+1, GetName());
247  }
248  }
249  if (fZ[isec]>fZ[isec+1]) {
250  InspectShape();
251  Fatal("ComputeBBox", "Wrong section order");
252  }
253  }
254  // Check if the last sections are valid
255  if (TMath::Abs(fZ[1]-fZ[0]) < TGeoShape::Tolerance() ||
257  InspectShape();
258  Fatal("ComputeBBox","Shape %s at index %d: Not allowed first two or last two sections at same Z",
260  }
261  Double_t zmin = TMath::Min(fZ[0], fZ[fNz-1]);
262  Double_t zmax = TMath::Max(fZ[0], fZ[fNz-1]);
263  // find largest rmax an smallest rmin
264  Double_t rmin, rmax;
265  rmin = fRmin[TMath::LocMin(fNz, fRmin)];
266  rmax = fRmax[TMath::LocMax(fNz, fRmax)];
267 
268  Double_t xc[4];
269  Double_t yc[4];
270  xc[0] = rmax*fC1;
271  yc[0] = rmax*fS1;
272  xc[1] = rmax*fC2;
273  yc[1] = rmax*fS2;
274  xc[2] = rmin*fC1;
275  yc[2] = rmin*fS1;
276  xc[3] = rmin*fC2;
277  yc[3] = rmin*fS2;
278 
279  Double_t xmin = xc[TMath::LocMin(4, &xc[0])];
280  Double_t xmax = xc[TMath::LocMax(4, &xc[0])];
281  Double_t ymin = yc[TMath::LocMin(4, &yc[0])];
282  Double_t ymax = yc[TMath::LocMax(4, &yc[0])];
283 
284  Double_t ddp = -fPhi1;
285  if (ddp<0) ddp+= 360;
286  if (ddp<=fDphi) xmax = rmax;
287  ddp = 90-fPhi1;
288  if (ddp<0) ddp+= 360;
289  if (ddp<=fDphi) ymax = rmax;
290  ddp = 180-fPhi1;
291  if (ddp<0) ddp+= 360;
292  if (ddp<=fDphi) xmin = -rmax;
293  ddp = 270-fPhi1;
294  if (ddp<0) ddp+= 360;
295  if (ddp<=fDphi) ymin = -rmax;
296  fOrigin[0] = (xmax+xmin)/2;
297  fOrigin[1] = (ymax+ymin)/2;
298  fOrigin[2] = (zmax+zmin)/2;
299  fDX = (xmax-xmin)/2;
300  fDY = (ymax-ymin)/2;
301  fDZ = (zmax-zmin)/2;
303 }
304 
305 ////////////////////////////////////////////////////////////////////////////////
306 /// Compute normal to closest surface from POINT.
307 
308 void TGeoPcon::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
309 {
310  memset(norm,0,3*sizeof(Double_t));
311  Double_t r;
312  Double_t ptnew[3];
313  Double_t dz, rmin1, rmax1, rmin2, rmax2;
314  Bool_t is_tube;
315  Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]);
316  if (ipl==(fNz-1) || ipl<0) {
317  // point outside Z range
318  norm[2] = TMath::Sign(1., dir[2]);
319  return;
320  }
321  Int_t iplclose = ipl;
322  if ((fZ[ipl+1]-point[2])<(point[2]-fZ[ipl])) iplclose++;
323  dz = TMath::Abs(fZ[iplclose]-point[2]);
324  if (dz<1E-5) {
325  if (iplclose==0 || iplclose==(fNz-1)) {
326  norm[2] = TMath::Sign(1., dir[2]);
327  return;
328  }
329  if (iplclose==ipl && TGeoShape::IsSameWithinTolerance(fZ[ipl],fZ[ipl-1])) {
330  r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
331  if (r<TMath::Max(fRmin[ipl],fRmin[ipl-1]) || r>TMath::Min(fRmax[ipl],fRmax[ipl-1])) {
332  norm[2] = TMath::Sign(1., dir[2]);
333  return;
334  }
335  } else {
336  if (TGeoShape::IsSameWithinTolerance(fZ[iplclose],fZ[iplclose+1])) {
337  r = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
338  if (r<TMath::Max(fRmin[iplclose],fRmin[iplclose+1]) || r>TMath::Min(fRmax[iplclose],fRmax[iplclose+1])) {
339  norm[2] = TMath::Sign(1., dir[2]);
340  return;
341  }
342  }
343  }
344  } //-> Z done
345  memcpy(ptnew, point, 3*sizeof(Double_t));
346  dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
348  norm[2] = TMath::Sign(1., dir[2]);
349  return;
350  }
351  ptnew[2] -= 0.5*(fZ[ipl]+fZ[ipl+1]);
352  rmin1 = fRmin[ipl];
353  rmax1 = fRmax[ipl];
354  rmin2 = fRmin[ipl+1];
355  rmax2 = fRmax[ipl+1];
356  is_tube = (TGeoShape::IsSameWithinTolerance(rmin1,rmin2) && TGeoShape::IsSameWithinTolerance(rmax1,rmax2))?kTRUE:kFALSE;
357  if (!fFullPhi) {
358  if (is_tube) TGeoTubeSeg::ComputeNormalS(ptnew,dir,norm,rmin1,rmax1,dz,fC1,fS1,fC2,fS2);
359  else TGeoConeSeg::ComputeNormalS(ptnew,dir,norm,dz,rmin1,rmax1,rmin2,rmax2,fC1,fS1,fC2,fS2);
360  } else {
361  if (is_tube) TGeoTube::ComputeNormalS(ptnew,dir,norm,rmin1,rmax1,dz);
362  else TGeoCone::ComputeNormalS(ptnew,dir,norm,dz,rmin1,rmax1,rmin2,rmax2);
363  }
364 }
365 
366 ////////////////////////////////////////////////////////////////////////////////
367 /// test if point is inside this shape
368 /// check total z range
369 
370 Bool_t TGeoPcon::Contains(const Double_t *point) const
371 {
372  if ((point[2]<fZ[0]) || (point[2]>fZ[fNz-1])) return kFALSE;
373  // check R squared
374  Double_t r2 = point[0]*point[0]+point[1]*point[1];
375 
376  Int_t izl = 0;
377  Int_t izh = fNz-1;
378  Int_t izt = (fNz-1)/2;
379  while ((izh-izl)>1) {
380  if (point[2] > fZ[izt]) izl = izt;
381  else izh = izt;
382  izt = (izl+izh)>>1;
383  }
384  // the point is in the section bounded by izl and izh Z planes
385 
386  // compute Rmin and Rmax and test the value of R squared
387  Double_t rmin, rmax;
388  if (TGeoShape::IsSameWithinTolerance(fZ[izl],fZ[izh]) && TGeoShape::IsSameWithinTolerance(point[2],fZ[izl])) {
389  rmin = TMath::Min(fRmin[izl], fRmin[izh]);
390  rmax = TMath::Max(fRmax[izl], fRmax[izh]);
391  } else {
392  Double_t dz = fZ[izh] - fZ[izl];
393  Double_t dz1 = point[2] - fZ[izl];
394  rmin = (fRmin[izl]*(dz-dz1)+fRmin[izh]*dz1)/dz;
395  rmax = (fRmax[izl]*(dz-dz1)+fRmax[izh]*dz1)/dz;
396  }
397  if ((r2<rmin*rmin) || (r2>rmax*rmax)) return kFALSE;
398  // now check phi
400  if (r2<1E-10) return kTRUE;
401  Double_t phi = TMath::ATan2(point[1], point[0]) * TMath::RadToDeg();
402  if (phi < 0) phi+=360.0;
403  Double_t ddp = phi-fPhi1;
404  if (ddp<0) ddp+=360.;
405  if (ddp<=fDphi) return kTRUE;
406  return kFALSE;
407 }
408 
409 ////////////////////////////////////////////////////////////////////////////////
410 /// compute closest distance from point px,py to each corner
411 
413 {
415  const Int_t numPoints = 2*n*fNz;
416  return ShapeDistancetoPrimitive(numPoints, px, py);
417 }
418 
419 ////////////////////////////////////////////////////////////////////////////////
420 /// compute distance from inside point to surface of the polycone
421 
422 Double_t TGeoPcon::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
423 {
424  if (iact<3 && safe) {
425  *safe = Safety(point, kTRUE);
426  if (iact==0) return TGeoShape::Big();
427  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
428  }
429  Double_t snxt = TGeoShape::Big();
430  Double_t sstep = 1E-6;
431  Double_t point_new[3];
432  // determine which z segment contains the point
433  Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]+TMath::Sign(1.E-10,dir[2]));
434  if (ipl<0) ipl=0;
435  if (ipl==(fNz-1)) ipl--;
436  Double_t dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
437  Bool_t special_case = kFALSE;
438  if (dz<1e-9) {
439  // radius changing segment, make sure track is not in the XY plane
440  if (TGeoShape::IsSameWithinTolerance(dir[2], 0)) {
441  special_case = kTRUE;
442  } else {
443  //check if a close point is still contained
444  point_new[0] = point[0]+sstep*dir[0];
445  point_new[1] = point[1]+sstep*dir[1];
446  point_new[2] = point[2]+sstep*dir[2];
447  if (!Contains(point_new)) return 0.;
448  return (DistFromInside(point_new,dir,iact,step,safe)+sstep);
449  }
450  }
451  // determine if the current segment is a tube or a cone
452  Bool_t intub = kTRUE;
453  if (!TGeoShape::IsSameWithinTolerance(fRmin[ipl],fRmin[ipl+1])) intub=kFALSE;
454  else if (!TGeoShape::IsSameWithinTolerance(fRmax[ipl],fRmax[ipl+1])) intub=kFALSE;
455  // determine phi segmentation
456  memcpy(point_new, point, 2*sizeof(Double_t));
457  // new point in reference system of the current segment
458  point_new[2] = point[2]-0.5*(fZ[ipl]+fZ[ipl+1]);
459 
460  if (special_case) {
461  if (!fFullPhi)
462  snxt = TGeoTubeSeg::DistFromInsideS(point_new, dir,
463  TMath::Min(fRmin[ipl],fRmin[ipl+1]), TMath::Max(fRmax[ipl],fRmax[ipl+1]),
464  dz, fC1,fS1,fC2,fS2,fCm,fSm,fCdphi);
465  else
466  snxt = TGeoTube::DistFromInsideS(point_new, dir,
467  TMath::Min(fRmin[ipl],fRmin[ipl+1]), TMath::Max(fRmax[ipl],fRmax[ipl+1]),dz);
468  return snxt;
469  }
470  if (intub) {
471  if (!fFullPhi)
472  snxt = TGeoTubeSeg::DistFromInsideS(point_new, dir, fRmin[ipl], fRmax[ipl],dz, fC1,fS1,fC2,fS2,fCm,fSm,fCdphi);
473  else
474  snxt = TGeoTube::DistFromInsideS(point_new, dir, fRmin[ipl], fRmax[ipl],dz);
475  } else {
476  if (!fFullPhi)
477  snxt = TGeoConeSeg::DistFromInsideS(point_new,dir,dz,fRmin[ipl],fRmax[ipl],fRmin[ipl+1],fRmax[ipl+1],fC1,fS1,fC2,fS2,fCm,fSm,fCdphi);
478  else
479  snxt=TGeoCone::DistFromInsideS(point_new,dir,dz,fRmin[ipl],fRmax[ipl],fRmin[ipl+1], fRmax[ipl+1]);
480  }
481 
482  for (Int_t i=0; i<3; i++) point_new[i]=point[i]+(snxt+1E-6)*dir[i];
483  if (!Contains(&point_new[0])) return snxt;
484 
485  snxt += DistFromInside(&point_new[0], dir, 3) + 1E-6;
486  return snxt;
487 }
488 
489 ////////////////////////////////////////////////////////////////////////////////
490 /// compute distance to a pcon Z slice. Segment iz must be valid
491 
492 Double_t TGeoPcon::DistToSegZ(const Double_t *point, const Double_t *dir, Int_t &iz) const
493 {
494  Double_t zmin=fZ[iz];
495  Double_t zmax=fZ[iz+1];
496  if (TGeoShape::IsSameWithinTolerance(zmin,zmax)) {
497  if (TGeoShape::IsSameWithinTolerance(dir[2],0)) return TGeoShape::Big();
498  Int_t istep=(dir[2]>0)?1:-1;
499  iz+=istep;
500  if (iz<0 || iz>(fNz-2)) return TGeoShape::Big();
501  return DistToSegZ(point,dir,iz);
502  }
503  Double_t dz=0.5*(zmax-zmin);
504  Double_t local[3];
505  memcpy(&local[0], point, 3*sizeof(Double_t));
506  local[2]=point[2]-0.5*(zmin+zmax);
507  Double_t snxt;
508  Double_t rmin1=fRmin[iz];
509  Double_t rmax1=fRmax[iz];
510  Double_t rmin2=fRmin[iz+1];
511  Double_t rmax2=fRmax[iz+1];
512 
513  if (TGeoShape::IsSameWithinTolerance(rmin1,rmin2) && TGeoShape::IsSameWithinTolerance(rmax1,rmax2)) {
514  if (fFullPhi) snxt=TGeoTube::DistFromOutsideS(local, dir, rmin1, rmax1, dz);
515  else snxt=TGeoTubeSeg::DistFromOutsideS(local,dir,rmin1,rmax1,dz,fC1,fS1,fC2,fS2,fCm,fSm,fCdphi);
516  } else {
517  if (fFullPhi) snxt=TGeoCone::DistFromOutsideS(local,dir,dz,rmin1, rmax1,rmin2,rmax2);
518  else snxt=TGeoConeSeg::DistFromOutsideS(local,dir,dz,rmin1,rmax1,rmin2,rmax2,fC1,fS1,fC2,fS2,fCm,fSm,fCdphi);
519  }
520  if (snxt<1E20) return snxt;
521  // check next segment
522  if (TGeoShape::IsSameWithinTolerance(dir[2],0)) return TGeoShape::Big();
523  Int_t istep=(dir[2]>0)?1:-1;
524  iz+=istep;
525  if (iz<0 || iz>(fNz-2)) return TGeoShape::Big();
526  return DistToSegZ(point,dir,iz);
527 }
528 
529 ////////////////////////////////////////////////////////////////////////////////
530 /// compute distance from outside point to surface of the tube
531 
532 Double_t TGeoPcon::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
533 {
534  if ((iact<3) && safe) {
535  *safe = Safety(point, kFALSE);
536  if ((iact==1) && (*safe>step)) return TGeoShape::Big();
537  if (iact==0) return TGeoShape::Big();
538  }
539  // check if ray intersect outscribed cylinder
540  if ((point[2]<fZ[0]) && (dir[2]<=0)) return TGeoShape::Big();
541  if ((point[2]>fZ[fNz-1]) && (dir[2]>=0)) return TGeoShape::Big();
542 // Check if the bounding box is crossed within the requested distance
543  Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
544  if (sdist>=step) return TGeoShape::Big();
545 
546  Double_t r2 = point[0]*point[0]+point[1]*point[1];
547  Double_t radmax=0;
548  radmax=fRmax[TMath::LocMax(fNz, fRmax)];
549  if (r2>(radmax*radmax)) {
550  Double_t rpr=-point[0]*dir[0]-point[1]*dir[1];
551  Double_t nxy=dir[0]*dir[0]+dir[1]*dir[1];
552  if (rpr<TMath::Sqrt((r2-radmax*radmax)*nxy)) return TGeoShape::Big();
553  }
554 
555  // find in which Z segment we are
556  Int_t ipl = TMath::BinarySearch(fNz, fZ, point[2]);
557  Int_t ifirst = ipl;
558  if (ifirst<0) {
559  ifirst = 0;
560  } else if (ifirst>=(fNz-1)) {
561  ifirst = fNz-2;
562  }
563  // find if point is in the phi gap
564  Double_t phi=0;
565  if (!fFullPhi) {
566  phi = TMath::ATan2(point[1], point[0]);
567  if (phi<0) phi+=2.*TMath::Pi();
568  }
569 
570  // compute distance to boundary
571  return DistToSegZ(point,dir,ifirst);
572 }
573 
574 ////////////////////////////////////////////////////////////////////////////////
575 /// Defines z position of a section plane, rmin and rmax at this z. Sections
576 /// should be defined in increasing or decreasing Z order and the last section
577 /// HAS to be snum = fNz-1
578 
580 {
581  if ((snum<0) || (snum>=fNz)) return;
582  fZ[snum] = z;
583  fRmin[snum] = rmin;
584  fRmax[snum] = rmax;
585  if (rmin>rmax)
586  Warning("DefineSection", "Shape %s: invalid rmin=%g rmax=%g", GetName(), rmin, rmax);
587  if (snum==(fNz-1)) {
588  // Reorder sections in increasing Z order
589  if (fZ[0] > fZ[snum]) {
590  Int_t iz = 0;
591  Int_t izi = fNz-1;
592  Double_t temp;
593  while (iz<izi) {
594  temp = fZ[iz];
595  fZ[iz] = fZ[izi];
596  fZ[izi] = temp;
597  temp = fRmin[iz];
598  fRmin[iz] = fRmin[izi];
599  fRmin[izi] = temp;
600  temp = fRmax[iz];
601  fRmax[iz] = fRmax[izi];
602  fRmax[izi] = temp;
603  iz++;
604  izi--;
605  }
606  }
607  ComputeBBox();
608  }
609 }
610 
611 ////////////////////////////////////////////////////////////////////////////////
612 /// Returns number of segments on each mesh circle segment.
613 
615 {
616  return gGeoManager->GetNsegments();
617 }
618 
619 ////////////////////////////////////////////////////////////////////////////////
620 /// Divide this polycone shape belonging to volume "voldiv" into ndiv volumes
621 /// called divname, from start position with the given step. Returns pointer
622 /// to created division cell volume in case of Z divisions. Z divisions can be
623 /// performed if the divided range is in between two consecutive Z planes.
624 /// In case a wrong division axis is supplied, returns pointer to
625 /// volume that was divided.
626 
627 TGeoVolume *TGeoPcon::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
628  Double_t start, Double_t step)
629 {
630  TGeoShape *shape; //--- shape to be created
631  TGeoVolume *vol; //--- division volume to be created
632  TGeoVolumeMulti *vmulti; //--- generic divided volume
633  TGeoPatternFinder *finder; //--- finder to be attached
634  TString opt = ""; //--- option to be attached
635  Double_t zmin = start;
636  Double_t zmax = start+ndiv*step;
637  Int_t isect = -1;
638  Int_t is, id, ipl;
639  switch (iaxis) {
640  case 1: //--- R division
641  Error("Divide", "Shape %s: cannot divide a pcon on radius", GetName());
642  return 0;
643  case 2: //--- Phi division
644  finder = new TGeoPatternCylPhi(voldiv, ndiv, start, start+ndiv*step);
645  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
646  voldiv->SetFinder(finder);
647  finder->SetDivIndex(voldiv->GetNdaughters());
648  shape = new TGeoPcon(-step/2, step, fNz);
649  for (is=0; is<fNz; is++)
650  ((TGeoPcon*)shape)->DefineSection(is, fZ[is], fRmin[is], fRmax[is]);
651  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
652  vmulti->AddVolume(vol);
653  opt = "Phi";
654  for (id=0; id<ndiv; id++) {
655  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
656  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
657  }
658  return vmulti;
659  case 3: //--- Z division
660  // find start plane
661  for (ipl=0; ipl<fNz-1; ipl++) {
662  if (start<fZ[ipl]) continue;
663  else {
664  if ((start+ndiv*step)>fZ[ipl+1]) continue;
665  }
666  isect = ipl;
667  zmin = fZ[isect];
668  zmax= fZ[isect+1];
669  break;
670  }
671  if (isect<0) {
672  Error("Divide", "Shape %s: cannot divide pcon on Z if divided region is not between 2 planes", GetName());
673  return 0;
674  }
675  finder = new TGeoPatternZ(voldiv, ndiv, start, start+ndiv*step);
676  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
677  voldiv->SetFinder(finder);
678  finder->SetDivIndex(voldiv->GetNdaughters());
679  opt = "Z";
680  for (id=0; id<ndiv; id++) {
681  Double_t z1 = start+id*step;
682  Double_t z2 = start+(id+1)*step;
683  Double_t rmin1 = (fRmin[isect]*(zmax-z1)-fRmin[isect+1]*(zmin-z1))/(zmax-zmin);
684  Double_t rmax1 = (fRmax[isect]*(zmax-z1)-fRmax[isect+1]*(zmin-z1))/(zmax-zmin);
685  Double_t rmin2 = (fRmin[isect]*(zmax-z2)-fRmin[isect+1]*(zmin-z2))/(zmax-zmin);
686  Double_t rmax2 = (fRmax[isect]*(zmax-z2)-fRmax[isect+1]*(zmin-z2))/(zmax-zmin);
688  Bool_t is_seg = (fDphi<360)?kTRUE:kFALSE;
689  if (is_seg) {
690  if (is_tube) shape=new TGeoTubeSeg(fRmin[isect],fRmax[isect],step/2, fPhi1, fPhi1+fDphi);
691  else shape=new TGeoConeSeg(step/2, rmin1, rmax1, rmin2, rmax2, fPhi1, fPhi1+fDphi);
692  } else {
693  if (is_tube) shape=new TGeoTube(fRmin[isect],fRmax[isect],step/2);
694  else shape = new TGeoCone(step/2,rmin1,rmax1,rmin2,rmax2);
695  }
696  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
697  vmulti->AddVolume(vol);
698  voldiv->AddNodeOffset(vol, id, start+id*step+step/2, opt.Data());
699  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
700  }
701  return vmulti;
702  default:
703  Error("Divide", "Shape %s: Wrong axis %d for division",GetName(), iaxis);
704  return 0;
705  }
706 }
707 
708 ////////////////////////////////////////////////////////////////////////////////
709 /// Returns name of axis IAXIS.
710 
711 const char *TGeoPcon::GetAxisName(Int_t iaxis) const
712 {
713  switch (iaxis) {
714  case 1:
715  return "R";
716  case 2:
717  return "PHI";
718  case 3:
719  return "Z";
720  default:
721  return "UNDEFINED";
722  }
723 }
724 
725 ////////////////////////////////////////////////////////////////////////////////
726 /// Get range of shape for a given axis.
727 
729 {
730  xlo = 0;
731  xhi = 0;
732  Double_t dx = 0;
733  switch (iaxis) {
734  case 2:
735  xlo = fPhi1;
736  xhi = fPhi1 + fDphi;
737  dx = fDphi;
738  return dx;
739  case 3:
740  xlo = fZ[0];
741  xhi = fZ[fNz-1];
742  dx = xhi-xlo;
743  return dx;
744  }
745  return dx;
746 }
747 
748 ////////////////////////////////////////////////////////////////////////////////
749 /// Fill vector param[4] with the bounding cylinder parameters. The order
750 /// is the following : Rmin, Rmax, Phi1, Phi2
751 
753 {
754  param[0] = fRmin[0]; // Rmin
755  param[1] = fRmax[0]; // Rmax
756  for (Int_t i=1; i<fNz; i++) {
757  if (fRmin[i] < param[0]) param[0] = fRmin[i];
758  if (fRmax[i] > param[1]) param[1] = fRmax[i];
759  }
760  param[0] *= param[0];
761  param[1] *= param[1];
763  param[2] = 0.;
764  param[3] = 360.;
765  return;
766  }
767  param[2] = (fPhi1<0)?(fPhi1+360.):fPhi1; // Phi1
768  param[3] = param[2]+fDphi; // Phi2
769 }
770 
771 ////////////////////////////////////////////////////////////////////////////////
772 /// Returns Rmin for Z segment IPL.
773 
775 {
776  if (ipl<0 || ipl>(fNz-1)) {
777  Error("GetRmin","ipl=%i out of range (0,%i) in shape %s",ipl,fNz-1,GetName());
778  return 0.;
779  }
780  return fRmin[ipl];
781 }
782 
783 ////////////////////////////////////////////////////////////////////////////////
784 /// Returns Rmax for Z segment IPL.
785 
787 {
788  if (ipl<0 || ipl>(fNz-1)) {
789  Error("GetRmax","ipl=%i out of range (0,%i) in shape %s",ipl,fNz-1,GetName());
790  return 0.;
791  }
792  return fRmax[ipl];
793 }
794 
795 ////////////////////////////////////////////////////////////////////////////////
796 /// Returns Z for segment IPL.
797 
799 {
800  if (ipl<0 || ipl>(fNz-1)) {
801  Error("GetZ","ipl=%i out of range (0,%i) in shape %s",ipl,fNz-1,GetName());
802  return 0.;
803  }
804  return fZ[ipl];
805 }
806 
807 ////////////////////////////////////////////////////////////////////////////////
808 /// print shape parameters
809 
811 {
812  printf("*** Shape %s: TGeoPcon ***\n", GetName());
813  printf(" Nz = %i\n", fNz);
814  printf(" phi1 = %11.5f\n", fPhi1);
815  printf(" dphi = %11.5f\n", fDphi);
816  for (Int_t ipl=0; ipl<fNz; ipl++)
817  printf(" plane %i: z=%11.5f Rmin=%11.5f Rmax=%11.5f\n", ipl, fZ[ipl], fRmin[ipl], fRmax[ipl]);
818  printf(" Bounding box:\n");
820 }
821 
822 ////////////////////////////////////////////////////////////////////////////////
823 /// Creates a TBuffer3D describing *this* shape.
824 /// Coordinates are in local reference frame.
825 
827 {
828  Int_t nbPnts, nbSegs, nbPols;
829  GetMeshNumbers(nbPnts, nbSegs, nbPols);
830  if (nbPnts <= 0) return nullptr;
831 
833  nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols);
834  if (buff)
835  {
836  SetPoints(buff->fPnts);
837  SetSegsAndPols(*buff);
838  }
839 
840  return buff;
841 }
842 
843 ////////////////////////////////////////////////////////////////////////////////
844 /// Fill TBuffer3D structure for segments and polygons.
845 
847 {
848  if (!HasInsideSurface()) {
850  return;
851  }
852 
853  Int_t i, j;
854  const Int_t n = gGeoManager->GetNsegments()+1;
855  Int_t nz = GetNz();
856  if (nz < 2) return;
857  Int_t nbPnts = nz*2*n;
858  if (nbPnts <= 0) return;
859  Double_t dphi = GetDphi();
860 
861  Bool_t specialCase = TGeoShape::IsSameWithinTolerance(dphi,360);
862  Int_t c = GetBasicColor();
863 
864  Int_t indx = 0, indx2, k;
865 
866  //inside & outside circles, number of segments: 2*nz*(n-1)
867  // special case number of segments: 2*nz*n
868  for (i = 0; i < nz*2; i++) {
869  indx2 = i*n;
870  for (j = 1; j < n; j++) {
871  buff.fSegs[indx++] = c;
872  buff.fSegs[indx++] = indx2+j-1;
873  buff.fSegs[indx++] = indx2+j;
874  }
875  if (specialCase) {
876  buff.fSegs[indx++] = c;
877  buff.fSegs[indx++] = indx2+j-1;
878  buff.fSegs[indx++] = indx2;
879  }
880  }
881 
882  //bottom & top lines, number of segments: 2*n
883  for (i = 0; i < 2; i++) {
884  indx2 = i*(nz-1)*2*n;
885  for (j = 0; j < n; j++) {
886  buff.fSegs[indx++] = c;
887  buff.fSegs[indx++] = indx2+j;
888  buff.fSegs[indx++] = indx2+n+j;
889  }
890  }
891 
892  //inside & outside cylinders, number of segments: 2*(nz-1)*n
893  for (i = 0; i < (nz-1); i++) {
894  //inside cylinder
895  indx2 = i*n*2;
896  for (j = 0; j < n; j++) {
897  buff.fSegs[indx++] = c+2;
898  buff.fSegs[indx++] = indx2+j;
899  buff.fSegs[indx++] = indx2+n*2+j;
900  }
901  //outside cylinder
902  indx2 = i*n*2+n;
903  for (j = 0; j < n; j++) {
904  buff.fSegs[indx++] = c+3;
905  buff.fSegs[indx++] = indx2+j;
906  buff.fSegs[indx++] = indx2+n*2+j;
907  }
908  }
909 
910  //left & right sections, number of segments: 2*(nz-2)
911  // special case number of segments: 0
912  if (!specialCase) {
913  for (i = 1; i < (nz-1); i++) {
914  for (j = 0; j < 2; j++) {
915  buff.fSegs[indx++] = c;
916  buff.fSegs[indx++] = 2*i * n + j*(n-1);
917  buff.fSegs[indx++] = (2*i+1) * n + j*(n-1);
918  }
919  }
920  }
921 
922  Int_t m = n - 1 + (specialCase ? 1 : 0);
923  indx = 0;
924 
925  //bottom & top, number of polygons: 2*(n-1)
926  // special case number of polygons: 2*n
927  for (j = 0; j < n-1; j++) {
928  buff.fPols[indx++] = c+3;
929  buff.fPols[indx++] = 4;
930  buff.fPols[indx++] = 2*nz*m+j;
931  buff.fPols[indx++] = m+j;
932  buff.fPols[indx++] = 2*nz*m+j+1;
933  buff.fPols[indx++] = j;
934  }
935  for (j = 0; j < n-1; j++) {
936  buff.fPols[indx++] = c+3;
937  buff.fPols[indx++] = 4;
938  buff.fPols[indx++] = 2*nz*m+n+j;
939  buff.fPols[indx++] = (nz*2-2)*m+j;
940  buff.fPols[indx++] = 2*nz*m+n+j+1;
941  buff.fPols[indx++] = (nz*2-2)*m+m+j;
942  }
943  if (specialCase) {
944  buff.fPols[indx++] = c+3;
945  buff.fPols[indx++] = 4;
946  buff.fPols[indx++] = 2*nz*m+j;
947  buff.fPols[indx++] = m+j;
948  buff.fPols[indx++] = 2*nz*m;
949  buff.fPols[indx++] = j;
950 
951  buff.fPols[indx++] = c+3;
952  buff.fPols[indx++] = 4;
953  buff.fPols[indx++] = 2*nz*m+n+j;
954  buff.fPols[indx++] = (nz*2-2)*m+m+j;
955  buff.fPols[indx++] = 2*nz*m+n;
956  buff.fPols[indx++] = (nz*2-2)*m+j;
957  }
958 
959  //inside & outside, number of polygons: (nz-1)*2*(n-1)
960  for (k = 0; k < (nz-1); k++) {
961  for (j = 0; j < n-1; j++) {
962  buff.fPols[indx++] = c;
963  buff.fPols[indx++] = 4;
964  buff.fPols[indx++] = 2*k*m+j;
965  buff.fPols[indx++] = nz*2*m+(2*k+2)*n+j+1;
966  buff.fPols[indx++] = (2*k+2)*m+j;
967  buff.fPols[indx++] = nz*2*m+(2*k+2)*n+j;
968  }
969  for (j = 0; j < n-1; j++) {
970  buff.fPols[indx++] = c+1;
971  buff.fPols[indx++] = 4;
972  buff.fPols[indx++] = (2*k+1)*m+j;
973  buff.fPols[indx++] = nz*2*m+(2*k+3)*n+j;
974  buff.fPols[indx++] = (2*k+3)*m+j;
975  buff.fPols[indx++] = nz*2*m+(2*k+3)*n+j+1;
976  }
977  if (specialCase) {
978  buff.fPols[indx++] = c;
979  buff.fPols[indx++] = 4;
980  buff.fPols[indx++] = 2*k*m+j;
981  buff.fPols[indx++] = nz*2*m+(2*k+2)*n;
982  buff.fPols[indx++] = (2*k+2)*m+j;
983  buff.fPols[indx++] = nz*2*m+(2*k+2)*n+j;
984 
985  buff.fPols[indx++] = c+1;
986  buff.fPols[indx++] = 4;
987  buff.fPols[indx++] = (2*k+1)*m+j;
988  buff.fPols[indx++] = nz*2*m+(2*k+3)*n+j;
989  buff.fPols[indx++] = (2*k+3)*m+j;
990  buff.fPols[indx++] = nz*2*m+(2*k+3)*n;
991  }
992  }
993 
994  //left & right sections, number of polygons: 2*(nz-1)
995  // special case number of polygons: 0
996  if (!specialCase) {
997  indx2 = nz*2*(n-1);
998  for (k = 0; k < (nz-1); k++) {
999  buff.fPols[indx++] = c+2;
1000  buff.fPols[indx++] = 4;
1001  buff.fPols[indx++] = k==0 ? indx2 : indx2+2*nz*n+2*(k-1);
1002  buff.fPols[indx++] = indx2+2*(k+1)*n;
1003  buff.fPols[indx++] = indx2+2*nz*n+2*k;
1004  buff.fPols[indx++] = indx2+(2*k+3)*n;
1005 
1006  buff.fPols[indx++] = c+2;
1007  buff.fPols[indx++] = 4;
1008  buff.fPols[indx++] = k==0 ? indx2+n-1 : indx2+2*nz*n+2*(k-1)+1;
1009  buff.fPols[indx++] = indx2+(2*k+3)*n+n-1;
1010  buff.fPols[indx++] = indx2+2*nz*n+2*k+1;
1011  buff.fPols[indx++] = indx2+2*(k+1)*n+n-1;
1012  }
1013  buff.fPols[indx-8] = indx2+n;
1014  buff.fPols[indx-2] = indx2+2*n-1;
1015  }
1016 }
1017 
1018 
1019 ////////////////////////////////////////////////////////////////////////////////
1020 /// Fill TBuffer3D structure for segments and polygons, when no inner surface exists
1021 
1023 {
1024  const Int_t n = gGeoManager->GetNsegments() + 1;
1025  const Int_t nz = GetNz();
1026  const Int_t nbPnts = nz * n + 2;
1027 
1028  if ((nz < 2) || (nbPnts <= 0) || (n < 2)) return;
1029 
1030  Int_t c = GetBasicColor();
1031 
1032  Int_t indx = 0, indx1 = 0, indx2 = 0, i, j;
1033 
1034  // outside circles, number of segments: nz*n
1035  for (i = 0; i < nz; i++) {
1036  indx2 = i * n;
1037  for (j = 1; j < n; j++) {
1038  buff.fSegs[indx++] = c;
1039  buff.fSegs[indx++] = indx2 + j - 1;
1040  buff.fSegs[indx++] = indx2 + j % (n-1);
1041  }
1042  }
1043 
1044  indx2 = 0;
1045  // bottom lines
1046  for (j = 0; j < n; j++) {
1047  buff.fSegs[indx++] = c;
1048  buff.fSegs[indx++] = indx2 + j % (n-1);
1049  buff.fSegs[indx++] = nbPnts - 2;
1050  }
1051 
1052  indx2 = (nz-1) * n;
1053  // top lines
1054  for (j = 0; j < n; j++) {
1055  buff.fSegs[indx++] = c;
1056  buff.fSegs[indx++] = indx2 + j % (n-1);
1057  buff.fSegs[indx++] = nbPnts - 1;
1058  }
1059 
1060  // outside cylinders, number of segments: (nz-1)*n
1061  for (i = 0; i < (nz - 1); i++) {
1062  // outside cylinder
1063  indx2 = i * n;
1064  for (j = 0; j < n; j++) {
1065  buff.fSegs[indx++] = c;
1066  buff.fSegs[indx++] = indx2 + j % (n-1);
1067  buff.fSegs[indx++] = indx2 + n + j % (n-1);
1068  }
1069  }
1070 
1071  indx = 0;
1072 
1073  // bottom cap
1074  indx1 = 0; // start of first z layer
1075  indx2 = nz*(n-1);
1076  for (j = 0; j < n - 1; j++) {
1077  buff.fPols[indx++] = c;
1078  buff.fPols[indx++] = 3;
1079  buff.fPols[indx++] = indx1 + j;
1080  buff.fPols[indx++] = indx2 + j + 1;
1081  buff.fPols[indx++] = indx2 + j;
1082  }
1083 
1084  // top cap
1085  indx1 = (nz-1)*(n-1); // start last z layer
1086  indx2 = nz*(n-1) + n;
1087  for (j = 0; j < n - 1; j++) {
1088  buff.fPols[indx++] = c;
1089  buff.fPols[indx++] = 3;
1090  buff.fPols[indx++] = indx1 + j; // last z layer
1091  buff.fPols[indx++] = indx2 + j;
1092  buff.fPols[indx++] = indx2 + j + 1;
1093  }
1094 
1095  // outside, number of polygons: (nz-1)*(n-1)
1096  for (Int_t k = 0; k < (nz - 1); k++) {
1097  indx1 = k*(n-1);
1098  indx2 = nz*(n-1) + n*2 + k*n;
1099  for (j = 0; j < n-1; j++) {
1100  buff.fPols[indx++] = c;
1101  buff.fPols[indx++] = 4;
1102  buff.fPols[indx++] = indx1 + j;
1103  buff.fPols[indx++] = indx2 + j;
1104  buff.fPols[indx++] = indx1 + j + (n-1);
1105  buff.fPols[indx++] = indx2 + j + 1;
1106  }
1107  }
1108 }
1109 
1110 ////////////////////////////////////////////////////////////////////////////////
1111 /// Compute safety from POINT to segment between planes ipl, ipl+1 within safmin.
1112 
1113 Double_t TGeoPcon::SafetyToSegment(const Double_t *point, Int_t ipl, Bool_t in, Double_t safmin) const
1114 {
1115  if (ipl<0 || ipl>fNz-2) return (safmin+1.); // error in input plane
1116 // Get info about segment.
1117  Double_t dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
1118  if (dz<1E-9) return 1E9; // radius-changing segment
1119  Double_t ptnew[3];
1120  memcpy(ptnew, point, 3*sizeof(Double_t));
1121  ptnew[2] -= 0.5*(fZ[ipl]+fZ[ipl+1]);
1122  Double_t safe = TMath::Abs(ptnew[2])-dz;
1123  if (safe>safmin) return TGeoShape::Big(); // means: stop checking further segments
1124  Double_t rmin1 = fRmin[ipl];
1125  Double_t rmax1 = fRmax[ipl];
1126  Double_t rmin2 = fRmin[ipl+1];
1127  Double_t rmax2 = fRmax[ipl+1];
1128  Bool_t is_tube = (TGeoShape::IsSameWithinTolerance(rmin1,rmin2) && TGeoShape::IsSameWithinTolerance(rmax1,rmax2))?kTRUE:kFALSE;
1129  if (!fFullPhi) {
1130  if (is_tube) safe = TGeoTubeSeg::SafetyS(ptnew,in,rmin1,rmax1, dz,fPhi1,fPhi1+fDphi,0);
1131  else safe = TGeoConeSeg::SafetyS(ptnew,in,dz,rmin1,rmax1,rmin2,rmax2,fPhi1,fPhi1+fDphi,0);
1132  } else {
1133  if (is_tube) safe = TGeoTube::SafetyS(ptnew,in,rmin1,rmax1,dz,0);
1134  else safe = TGeoCone::SafetyS(ptnew,in,dz,rmin1,rmax1,rmin2,rmax2,0);
1135  }
1136  if (safe<0) safe=0;
1137  return safe;
1138 }
1139 
1140 ////////////////////////////////////////////////////////////////////////////////
1141 /// computes the closest distance from given point to this shape, according
1142 /// to option. The matching point on the shape is stored in spoint.
1143 /// localize the Z segment
1144 
1145 Double_t TGeoPcon::Safety(const Double_t *point, Bool_t in) const
1146 {
1147  Double_t safmin, saftmp;
1148  Double_t dz;
1149  Int_t ipl, iplane;
1150 
1151  if (in) {
1152  //---> point is inside pcon
1153  ipl = TMath::BinarySearch(fNz, fZ, point[2]);
1154  if (ipl==(fNz-1)) return 0; // point on last Z boundary
1155  if (ipl<0) return 0; // point on first Z boundary
1156  if (ipl>0 && TGeoShape::IsSameWithinTolerance(fZ[ipl-1],fZ[ipl]) && TGeoShape::IsSameWithinTolerance(point[2],fZ[ipl-1])) ipl--;
1157  dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
1158  if (dz<1E-8) {
1159  // Point on a segment-changing plane
1160  safmin = TMath::Min(point[2]-fZ[ipl-1],fZ[ipl+2]-point[2]);
1161  saftmp = TGeoShape::Big();
1162  if (fDphi<360) saftmp = TGeoShape::SafetyPhi(point,in,fPhi1,fPhi1+fDphi);
1163  if (saftmp<safmin) safmin = saftmp;
1164  Double_t radius = TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
1165  if (fRmin[ipl]>0) safmin = TMath::Min(safmin, radius-fRmin[ipl]);
1166  if (fRmin[ipl+1]>0) safmin = TMath::Min(safmin, radius-fRmin[ipl+1]);
1167  safmin = TMath::Min(safmin, fRmax[ipl]-radius);
1168  safmin = TMath::Min(safmin, fRmax[ipl+1]-radius);
1169  if (safmin<0) safmin = 0;
1170  return safmin;
1171  }
1172  // Check safety for current segment
1173  safmin = SafetyToSegment(point, ipl);
1174  if (safmin>1E10) {
1175  // something went wrong - point is not inside current segment
1176  return 0.;
1177  }
1178  if (safmin<1E-6) return TMath::Abs(safmin); // point on radius-changing plane
1179  // check increasing iplanes
1180 /*
1181  iplane = ipl+1;
1182  saftmp = 0.;
1183  while ((iplane<fNz-1) && saftmp<1E10) {
1184  saftmp = TMath::Abs(SafetyToSegment(point,iplane,kFALSE,safmin));
1185  if (saftmp<safmin) safmin=saftmp;
1186  iplane++;
1187  }
1188  // now decreasing nplanes
1189  iplane = ipl-1;
1190  saftmp = 0.;
1191  while ((iplane>=0) && saftmp<1E10) {
1192  saftmp = TMath::Abs(SafetyToSegment(point,iplane,kFALSE,safmin));
1193  if (saftmp<safmin) safmin=saftmp;
1194  iplane--;
1195  }
1196 */
1197  return safmin;
1198  }
1199  //---> point is outside pcon
1200  ipl = TMath::BinarySearch(fNz, fZ, point[2]);
1201  if (ipl<0) ipl=0;
1202  else if (ipl==fNz-1) ipl=fNz-2;
1203  dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
1204  if (dz<1E-8 && (ipl+2<fNz)) {
1205  ipl++;
1206  dz = 0.5*(fZ[ipl+1]-fZ[ipl]);
1207  }
1208  // Check safety for current segment
1209  safmin = SafetyToSegment(point, ipl, kFALSE);
1210  if (safmin<1E-6) return TMath::Abs(safmin); // point on radius-changing plane
1211  saftmp = 0.;
1212  // check increasing iplanes
1213  iplane = ipl+1;
1214  saftmp = 0.;
1215  while ((iplane<fNz-1) && saftmp<1E10) {
1216  saftmp = TMath::Abs(SafetyToSegment(point,iplane,kFALSE,safmin));
1217  if (saftmp<safmin) safmin=saftmp;
1218  iplane++;
1219  }
1220  // now decreasing nplanes
1221  iplane = ipl-1;
1222  saftmp = 0.;
1223  while ((iplane>=0) && saftmp<1E10) {
1224  saftmp = TMath::Abs(SafetyToSegment(point,iplane,kFALSE,safmin));
1225  if (saftmp<safmin) safmin=saftmp;
1226  iplane--;
1227  }
1228  return safmin;
1229 }
1230 
1231 ////////////////////////////////////////////////////////////////////////////////
1232 /// Save a primitive as a C++ statement(s) on output stream "out".
1233 
1234 void TGeoPcon::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1235 {
1236  if (TObject::TestBit(kGeoSavePrimitive)) return;
1237  out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
1238  out << " phi1 = " << fPhi1 << ";" << std::endl;
1239  out << " dphi = " << fDphi << ";" << std::endl;
1240  out << " nz = " << fNz << ";" << std::endl;
1241  out << " TGeoPcon *pcon = new TGeoPcon(\"" << GetName() << "\",phi1,dphi,nz);" << std::endl;
1242  for (Int_t i=0; i<fNz; i++) {
1243  out << " z = " << fZ[i] << ";" << std::endl;
1244  out << " rmin = " << fRmin[i] << ";" << std::endl;
1245  out << " rmax = " << fRmax[i] << ";" << std::endl;
1246  out << " pcon->DefineSection(" << i << ", z,rmin,rmax);" << std::endl;
1247  }
1248  out << " TGeoShape *" << GetPointerName() << " = pcon;" << std::endl;
1250 }
1251 
1252 ////////////////////////////////////////////////////////////////////////////////
1253 /// Set polycone dimensions starting from an array.
1254 
1256 {
1257  fPhi1 = param[0];
1258  while (fPhi1<0) fPhi1 += 360.;
1259  fDphi = param[1];
1260  fNz = (Int_t)param[2];
1261  if (fNz<2) {
1262  Error("SetDimensions","Pcon %s: Number of Z sections must be > 2", GetName());
1263  return;
1264  }
1265  if (fRmin) delete [] fRmin;
1266  if (fRmax) delete [] fRmax;
1267  if (fZ) delete [] fZ;
1268  fRmin = new Double_t [fNz];
1269  fRmax = new Double_t [fNz];
1270  fZ = new Double_t [fNz];
1271  memset(fRmin, 0, fNz*sizeof(Double_t));
1272  memset(fRmax, 0, fNz*sizeof(Double_t));
1273  memset(fZ, 0, fNz*sizeof(Double_t));
1275  Double_t phi1 = fPhi1;
1276  Double_t phi2 = phi1+fDphi;
1277  Double_t phim = 0.5*(phi1+phi2);
1278  fC1 = TMath::Cos(phi1*TMath::DegToRad());
1279  fS1 = TMath::Sin(phi1*TMath::DegToRad());
1280  fC2 = TMath::Cos(phi2*TMath::DegToRad());
1281  fS2 = TMath::Sin(phi2*TMath::DegToRad());
1282  fCm = TMath::Cos(phim*TMath::DegToRad());
1283  fSm = TMath::Sin(phim*TMath::DegToRad());
1285 
1286  for (Int_t i=0; i<fNz; i++)
1287  DefineSection(i, param[3+3*i], param[4+3*i], param[5+3*i]);
1288 }
1289 
1290 ////////////////////////////////////////////////////////////////////////////////
1291 /// create polycone mesh points
1292 
1294 {
1295  Double_t phi, dphi;
1296  Int_t n = gGeoManager->GetNsegments() + 1;
1297  dphi = fDphi/(n-1);
1298  Int_t i, j;
1299  Int_t indx = 0;
1300 
1301  Bool_t hasInside = HasInsideSurface();
1302 
1303  if (points) {
1304  for (i = 0; i < fNz; i++) {
1305  if (hasInside)
1306  for (j = 0; j < n; j++) {
1307  phi = (fPhi1+j*dphi)*TMath::DegToRad();
1308  points[indx++] = fRmin[i] * TMath::Cos(phi);
1309  points[indx++] = fRmin[i] * TMath::Sin(phi);
1310  points[indx++] = fZ[i];
1311  }
1312  for (j = 0; j < n; j++) {
1313  phi = (fPhi1+j*dphi)*TMath::DegToRad();
1314  points[indx++] = fRmax[i] * TMath::Cos(phi);
1315  points[indx++] = fRmax[i] * TMath::Sin(phi);
1316  points[indx++] = fZ[i];
1317  }
1318  }
1319  if (!hasInside) {
1320  points[indx++] = 0;
1321  points[indx++] = 0;
1322  points[indx++] = fZ[0];
1323 
1324  points[indx++] = 0;
1325  points[indx++] = 0;
1326  points[indx++] = fZ[GetNz()-1];
1327  }
1328  }
1329 }
1330 
1331 ////////////////////////////////////////////////////////////////////////////////
1332 /// create polycone mesh points
1333 
1335 {
1336  Double_t phi, dphi;
1337  Int_t n = gGeoManager->GetNsegments() + 1;
1338  dphi = fDphi/(n-1);
1339  Int_t i, j;
1340  Int_t indx = 0;
1341 
1342  Bool_t hasInside = HasInsideSurface();
1343 
1344  if (points) {
1345  for (i = 0; i < fNz; i++) {
1346  if (hasInside)
1347  for (j = 0; j < n; j++) {
1348  phi = (fPhi1+j*dphi)*TMath::DegToRad();
1349  points[indx++] = fRmin[i] * TMath::Cos(phi);
1350  points[indx++] = fRmin[i] * TMath::Sin(phi);
1351  points[indx++] = fZ[i];
1352  }
1353  for (j = 0; j < n; j++) {
1354  phi = (fPhi1+j*dphi)*TMath::DegToRad();
1355  points[indx++] = fRmax[i] * TMath::Cos(phi);
1356  points[indx++] = fRmax[i] * TMath::Sin(phi);
1357  points[indx++] = fZ[i];
1358  }
1359  }
1360  if (!hasInside) {
1361  points[indx++] = 0;
1362  points[indx++] = 0;
1363  points[indx++] = fZ[0];
1364 
1365  points[indx++] = 0;
1366  points[indx++] = 0;
1367  points[indx++] = fZ[GetNz()-1];
1368  }
1369  }
1370 }
1371 ////////////////////////////////////////////////////////////////////////////////
1372 /// Return number of vertices of the mesh representation
1373 
1375 {
1376  Int_t nvert, nsegs, npols;
1377  GetMeshNumbers(nvert, nsegs, npols);
1378  return nvert;
1379 }
1380 
1381 ////////////////////////////////////////////////////////////////////////////////
1382 /// fill size of this 3-D object
1383 
1385 {
1386 }
1387 
1388 ////////////////////////////////////////////////////////////////////////////////
1389 /// Returns true when pgon has internal surface
1390 /// It will be only disabled when all Rmin values are 0
1391 
1393 {
1394  // only when full 360 is used, internal part can be excluded
1395  Bool_t specialCase = TGeoShape::IsSameWithinTolerance(GetDphi(), 360);
1396  if (!specialCase) return kTRUE;
1397 
1398  for (Int_t i = 0; i < GetNz(); i++)
1399  if (fRmin[i] > 0.)
1400  return kTRUE;
1401 
1402  return kFALSE;
1403 }
1404 
1405 
1406 ////////////////////////////////////////////////////////////////////////////////
1407 /// Returns numbers of vertices, segments and polygons composing the shape mesh.
1408 
1409 void TGeoPcon::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const
1410 {
1411  nvert = nsegs = npols = 0;
1412 
1414  Int_t nz = GetNz();
1415  if (nz < 2) return;
1416 
1417  if (HasInsideSurface()) {
1418  Bool_t specialCase = TGeoShape::IsSameWithinTolerance(GetDphi(),360);
1419  nvert = nz*2*n;
1420  nsegs = 4*(nz*n-1+(specialCase ? 1 : 0));
1421  npols = 2*(nz*n-1+(specialCase ? 1 : 0));
1422  } else {
1423  nvert = nz * n + 2;
1424  nsegs = nz * (n - 1) + n * 2 + (nz - 1) * n;
1425  npols = 2 * (n - 1) + (nz - 1) * (n - 1);
1426  }
1427 }
1428 
1429 ////////////////////////////////////////////////////////////////////////////////
1430 /// Fills a static 3D buffer and returns a reference.
1431 
1432 const TBuffer3D & TGeoPcon::GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
1433 {
1434  static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
1435 
1436  TGeoBBox::FillBuffer3D(buffer, reqSections, localFrame);
1437 
1438  if (reqSections & TBuffer3D::kRawSizes) {
1439  Int_t nbPnts, nbSegs, nbPols;
1440  GetMeshNumbers(nbPnts, nbSegs, nbPols);
1441  if (nbPnts > 0) {
1442  if (buffer.SetRawSizes(nbPnts, 3*nbPnts, nbSegs, 3*nbSegs, nbPols, 6*nbPols)) {
1444  }
1445  }
1446  }
1447  // TODO: Push down to TGeoShape?? Would have to do raw sizes set first..
1448  // can rest of TGeoShape be deferred until after this?
1449  if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
1450  SetPoints(buffer.fPnts);
1451  if (!buffer.fLocalFrame) {
1452  TransformPoints(buffer.fPnts, buffer.NbPnts());
1453  }
1454 
1455  SetSegsAndPols(buffer);
1457  }
1458 
1459  return buffer;
1460 }
1461 
1462 ////////////////////////////////////////////////////////////////////////////////
1463 /// Stream an object of class TGeoPcon.
1464 
1465 void TGeoPcon::Streamer(TBuffer &R__b)
1466 {
1467  if (R__b.IsReading()) {
1468  R__b.ReadClassBuffer(TGeoPcon::Class(),this);
1470  Double_t phi1 = fPhi1;
1471  Double_t phi2 = phi1+fDphi;
1472  Double_t phim = 0.5*(phi1+phi2);
1473  fC1 = TMath::Cos(phi1*TMath::DegToRad());
1474  fS1 = TMath::Sin(phi1*TMath::DegToRad());
1475  fC2 = TMath::Cos(phi2*TMath::DegToRad());
1476  fS2 = TMath::Sin(phi2*TMath::DegToRad());
1477  fCm = TMath::Cos(phim*TMath::DegToRad());
1478  fSm = TMath::Sin(phim*TMath::DegToRad());
1480  } else {
1481  R__b.WriteClassBuffer(TGeoPcon::Class(),this);
1482  }
1483 }
1484 
1485 ////////////////////////////////////////////////////////////////////////////////
1486 /// Check the inside status for each of the points in the array.
1487 /// Input: Array of point coordinates + vector size
1488 /// Output: Array of Booleans for the inside of each point
1489 
1490 void TGeoPcon::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
1491 {
1492  for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
1493 }
1494 
1495 ////////////////////////////////////////////////////////////////////////////////
1496 /// Compute the normal for an array o points so that norm.dot.dir is positive
1497 /// Input: Arrays of point coordinates and directions + vector size
1498 /// Output: Array of normal directions
1499 
1500 void TGeoPcon::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
1501 {
1502  for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
1503 }
1504 
1505 ////////////////////////////////////////////////////////////////////////////////
1506 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
1507 
1508 void TGeoPcon::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1509 {
1510  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
1511 }
1512 
1513 ////////////////////////////////////////////////////////////////////////////////
1514 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
1515 
1516 void TGeoPcon::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1517 {
1518  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
1519 }
1520 
1521 ////////////////////////////////////////////////////////////////////////////////
1522 /// Compute safe distance from each of the points in the input array.
1523 /// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1524 /// Output: Safety values
1525 
1526 void TGeoPcon::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
1527 {
1528  for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
1529 }
TGeoPcon::fS1
Double_t fS1
Cosine of phi1.
Definition: TGeoPcon.h:29
TGeoTube::SafetyS
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Int_t skipz=0)
computes the closest distance from given point to this shape, according to option.
Definition: TGeoTube.cxx:864
c
#define c(i)
Definition: RSha256.hxx:101
TGeoPcon.h
TGeoPcon::Contains
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this shape check total z range
Definition: TGeoPcon.cxx:370
TGeoPcon::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoPcon.cxx:1234
m
auto * m
Definition: textangle.C:8
n
const Int_t n
Definition: legend1.C:16
TGeoPcon::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: TGeoPcon.cxx:1145
ymax
float ymax
Definition: THbookFile.cxx:95
TGeoVolume::SetFinder
void SetFinder(TGeoPatternFinder *finder)
Definition: TGeoVolume.h:231
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
TGeoTube::DistFromInsideS
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz)
Compute distance from inside point to surface of the tube (static) Boundary safe algorithm.
Definition: TGeoTube.cxx:286
e
#define e(i)
Definition: RSha256.hxx:103
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
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:67
TMath::ATan2
Double_t ATan2(Double_t y, Double_t x)
Definition: TMath.h:679
TGeoPcon::DistToSegZ
Double_t DistToSegZ(const Double_t *point, const Double_t *dir, Int_t &iz) const
compute distance to a pcon Z slice. Segment iz must be valid
Definition: TGeoPcon.cxx:492
TMath::RadToDeg
constexpr Double_t RadToDeg()
Conversion from radian to degree:
Definition: TMath.h:73
TGeoVolume::GetNdaughters
Int_t GetNdaughters() const
Definition: TGeoVolume.h:349
Option_t
const char Option_t
Definition: RtypesCore.h:66
TGeoPcon::SafetyToSegment
Double_t SafetyToSegment(const Double_t *point, Int_t ipl, Bool_t in=kTRUE, Double_t safmin=TGeoShape::Big()) const
Compute safety from POINT to segment between planes ipl, ipl+1 within safmin.
Definition: TGeoPcon.cxx:1113
TBuffer3D::SetSectionsValid
void SetSectionsValid(UInt_t mask)
Definition: TBuffer3D.h:65
TGeoTubeSeg::ComputeNormalS
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2)
Compute normal to closest surface from POINT.
Definition: TGeoTube.cxx:1406
TMath::Max
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
TMath::Cos
Double_t Cos(Double_t)
Definition: TMath.h:643
TGeoPcon::fRmax
Double_t * fRmax
Definition: TGeoPcon.h:25
gGeoManager
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
TString::Data
const char * Data() const
Definition: TString.h:369
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
TGeoPcon::Sizeof3D
virtual void Sizeof3D() const
fill size of this 3-D object
Definition: TGeoPcon.cxx:1384
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
TGeoPatternCylPhi
Definition: TGeoPatternFinder.h:396
TBuffer::ReadClassBuffer
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
r
ROOT::R::TRInterface & r
Definition: Object.C:4
xmax
float xmax
Definition: THbookFile.cxx:95
TGeoPcon::~TGeoPcon
virtual ~TGeoPcon()
destructor
Definition: TGeoPcon.cxx:205
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
TGeoPcon::GetBoundingCylinder
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoPcon.cxx:752
TGeoPcon::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: TGeoPcon.cxx:1516
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
TGeoPcon::SetSegsAndPolsNoInside
void SetSegsAndPolsNoInside(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons, when no inner surface exists.
Definition: TGeoPcon.cxx:1022
TMath::DegToRad
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:81
TGeoBBox::fOrigin
Double_t fOrigin[3]
Definition: TGeoBBox.h:24
Float_t
float Float_t
Definition: RtypesCore.h:57
TObjArray::IndexOf
Int_t IndexOf(const TObject *obj) const
Definition: TObjArray.cxx:605
TGeoPcon::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: TGeoPcon.cxx:1500
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:921
Int_t
int Int_t
Definition: RtypesCore.h:45
TGeoVolume.h
TGeoPcon::GetZ
Double_t * GetZ() const
Definition: TGeoPcon.h:83
TGeoPcon::GetNmeshVertices
virtual Int_t GetNmeshVertices() const
Return number of vertices of the mesh representation.
Definition: TGeoPcon.cxx:1374
TGeoManager::GetNsegments
Int_t GetNsegments() const
Get number of segments approximating circles.
Definition: TGeoManager.cxx:3351
TBuffer3D::NbPnts
UInt_t NbPnts() const
Definition: TBuffer3D.h:80
TGeoNodeOffset
Node containing an offset.
Definition: TGeoNode.h:184
TGeoShape::SetShapeBit
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
Definition: TGeoShape.cxx:524
TGeoCone
Conical tube class.
Definition: TGeoCone.h:18
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
TGeoTubeSeg
A phi segment of a tube.
Definition: TGeoTube.h:92
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
TGeoShape::GetBasicColor
Int_t GetBasicColor() const
Get the basic color (0-7).
Definition: TGeoShape.cxx:673
TGeoPcon::GetBuffer3D
virtual const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
Definition: TGeoPcon.cxx:1432
TObjArray::At
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
TBuffer3D::fSegs
Int_t * fSegs
Definition: TBuffer3D.h:113
TString
Basic string class.
Definition: TString.h:136
TGeoPcon::GetAxisName
virtual const char * GetAxisName(Int_t iaxis) const
Returns name of axis IAXIS.
Definition: TGeoPcon.cxx:711
TMath::LocMax
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition: TMath.h:1000
TGeoPcon::DistancetoPrimitive
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
compute closest distance from point px,py to each corner
Definition: TGeoPcon.cxx:412
TGeoPcon::fPhi1
Double_t fPhi1
Definition: TGeoPcon.h:22
TGeoCone.h
TGeoTubeSeg::DistFromOutsideS
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz, Double_t c1, Double_t s1, Double_t c2, Double_t s2, Double_t cm, Double_t sm, Double_t cdfi)
Static method to compute distance to arbitrary tube segment from outside point Boundary safe algorith...
Definition: TGeoTube.cxx:1519
TGeoPcon::ComputeNormal
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoPcon.cxx:308
TGeoPatternFinder
Base finder class for patterns.
Definition: TGeoPatternFinder.h:32
bool
TGeoVolume::GetMedium
TGeoMedium * GetMedium() const
Definition: TGeoVolume.h:173
TGeoPcon::TGeoPcon
TGeoPcon()
dummy ctor
Definition: TGeoPcon.cxx:68
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:2420
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:360
TGeoPcon::ComputeBBox
virtual void ComputeBBox()
compute bounding box of the pcon Check if the sections are in increasing Z order
Definition: TGeoPcon.cxx:238
TGeoPcon::fC2
Double_t fC2
Sine of phi1.
Definition: TGeoPcon.h:30
TGeoPcon::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: TGeoPcon.cxx:1490
TMath::LocMin
Long64_t LocMin(Long64_t n, const T *a)
Return index of array with the minimum element.
Definition: TMath.h:972
TBuffer3D
Generic 3D primitive description class.
Definition: TBuffer3D.h:18
TMath::Pi
constexpr Double_t Pi()
Definition: TMath.h:37
TGeoPcon::GetAxisRange
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
Definition: TGeoPcon.cxx:728
TGeoShape
Base abstract class for all shapes.
Definition: TGeoShape.h:26
TBuffer.h
TGeoPcon::fCm
Double_t fCm
Sine of phi1+dphi.
Definition: TGeoPcon.h:32
TGeoPcon::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: TGeoPcon.cxx:1526
xmin
float xmin
Definition: THbookFile.cxx:95
TBuffer3DTypes::kGeneric
@ kGeneric
Definition: TBuffer3DTypes.h:24
TGeoBBox::InspectShape
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoBBox.cxx:790
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:22
TGeoVolume::GetNodes
TObjArray * GetNodes()
Definition: TGeoVolume.h:167
TGeoPcon::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: TGeoPcon.cxx:532
TGeoPcon::MakeBuffer3D
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
Definition: TGeoPcon.cxx:826
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
TGeoPcon::InspectShape
virtual void InspectShape() const
print shape parameters
Definition: TGeoPcon.cxx:810
TBuffer::WriteClassBuffer
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
TBuffer3D.h
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
TMath::Sign
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:165
TGeoShape::kGeoSavePrimitive
@ kGeoSavePrimitive
Definition: TGeoShape.h:65
TGeoPcon::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 polycone
Definition: TGeoPcon.cxx:422
TGeoPcon::GetNz
Int_t GetNz() const
Definition: TGeoPcon.h:77
TGeoPcon::GetRmin
Double_t * GetRmin() const
Definition: TGeoPcon.h:79
TGeoTube
Cylindrical tube class.
Definition: TGeoTube.h:18
TBuffer3DTypes.h
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:639
TBuffer3D::kRaw
@ kRaw
Definition: TBuffer3D.h:54
TGeoPatternFinder::SetDivIndex
void SetDivIndex(Int_t index)
Definition: TGeoPatternFinder.h:99
TGeoShape::GetName
virtual const char * GetName() const
Get the shape name.
Definition: TGeoShape.cxx:248
TGeoPcon
A polycone.
Definition: TGeoPcon.h:18
TGeoTubeSeg::SafetyS
static Double_t SafetyS(const Double_t *point, Bool_t in, Double_t rmin, Double_t rmax, Double_t dz, Double_t phi1, Double_t phi2, Int_t skipz=0)
Static method to compute the closest distance from given point to this shape.
Definition: TGeoTube.cxx:2098
TGeoBBox
Box class.
Definition: TGeoBBox.h:18
TGeoPcon::fDphi
Double_t fDphi
Definition: TGeoPcon.h:23
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:1030
TGeoPcon::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: TGeoPcon.cxx:1508
TGeoTubeSeg::DistFromInsideS
static Double_t DistFromInsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz, 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 (static) Boundary safe algorithm.
Definition: TGeoTube.cxx:1455
TObject::Warning
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:879
TGeoPcon::SetPoints
virtual void SetPoints(Double_t *points) const
create polycone mesh points
Definition: TGeoPcon.cxx:1293
TMath::Min
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
TGeoPcon::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: TGeoPcon.cxx:1409
TGeoPcon::GetRmax
Double_t * GetRmax() const
Definition: TGeoPcon.h:81
ymin
float ymin
Definition: THbookFile.cxx:95
TMath::BinarySearch
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMathBase.h:278
TGeoPcon::DefineSection
virtual void DefineSection(Int_t snum, Double_t z, Double_t rmin, Double_t rmax)
Defines z position of a section plane, rmin and rmax at this z.
Definition: TGeoPcon.cxx:579
TGeoManager.h
TGeoBBox::fDZ
Double_t fDZ
Definition: TGeoBBox.h:23
TGeoPcon::fRmin
Double_t * fRmin
Definition: TGeoPcon.h:24
TBuffer::IsReading
Bool_t IsReading() const
Definition: TBuffer.h:86
TGeoTube::ComputeNormalS
static void ComputeNormalS(const Double_t *point, const Double_t *dir, Double_t *norm, Double_t rmin, Double_t rmax, Double_t dz)
Compute normal to closest surface from POINT.
Definition: TGeoTube.cxx:245
TGeoShape::GetPointerName
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:699
TGeoShape::kGeoClosedShape
@ kGeoClosedShape
Definition: TGeoShape.h:60
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
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
TBuffer3D::fLocalFrame
Bool_t fLocalFrame
Definition: TBuffer3D.h:90
TGeoPcon::SetDimensions
virtual void SetDimensions(Double_t *param)
Set polycone dimensions starting from an array.
Definition: TGeoPcon.cxx:1255
TGeoManager::MakeVolumeMulti
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Definition: TGeoManager.cxx:3310
TGeoPcon::fFullPhi
Bool_t fFullPhi
Definition: TGeoPcon.h:27
TGeoPcon::fZ
Double_t * fZ
Definition: TGeoPcon.h:26
TGeoPcon::SetSegsAndPols
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fill TBuffer3D structure for segments and polygons.
Definition: TGeoPcon.cxx:846
TGeoVolumeMulti
Volume families.
Definition: TGeoVolume.h:254
TGeoTube::DistFromOutsideS
static Double_t DistFromOutsideS(const Double_t *point, const Double_t *dir, Double_t rmin, Double_t rmax, Double_t dz)
Static method to compute distance from outside point to a tube with given parameters Boundary safe al...
Definition: TGeoTube.cxx:346
TGeoPatternZ
Definition: TGeoPatternFinder.h:185
points
point * points
Definition: X3DBuffer.c:22
TGeoPcon::HasInsideSurface
Bool_t HasInsideSurface() const
Returns true when pgon has internal surface It will be only disabled when all Rmin values are 0.
Definition: TGeoPcon.cxx:1392
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
TGeoPcon::fS2
Double_t fS2
Cosine of phi1+dphi.
Definition: TGeoPcon.h:31
TGeoPcon::Capacity
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoPcon.cxx:215
TBuffer3D::kRawSizes
@ kRawSizes
Definition: TBuffer3D.h:53
name
char name[80]
Definition: TGX11.cxx:110
TGeoPcon::fNz
Int_t fNz
Definition: TGeoPcon.h:21
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
TBuffer3D::fPnts
Double_t * fPnts
Definition: TBuffer3D.h:112
TGeoConeSeg
A phi segment of a conical tube.
Definition: TGeoCone.h:102
TGeoPcon::GetDphi
Double_t GetDphi() const
Definition: TGeoPcon.h:76
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
TGeoPcon::GetNsegments
virtual Int_t GetNsegments() const
Returns number of segments on each mesh circle segment.
Definition: TGeoPcon.cxx:614
TGeoPcon::fCdphi
Double_t fCdphi
Sine of (phi1+phi2)/2.
Definition: TGeoPcon.h:34
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
TGeoTube.h
TGeoPcon::fC1
Double_t fC1
Full phi range flag.
Definition: TGeoPcon.h:28
Class
void Class()
Definition: Class.C:29
TGeoManager::GetListOfShapes
TObjArray * GetListOfShapes() const
Definition: TGeoManager.h:495
TGeoPcon::fSm
Double_t fSm
Cosine of (phi1+phi2)/2.
Definition: TGeoPcon.h:33
TBuffer3D::fPols
Int_t * fPols
Definition: TBuffer3D.h:114
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
TMath::E
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:96
TGeoVolume
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition: TGeoVolume.h:49
TGeoBBox::fDX
Double_t fDX
Definition: TGeoBBox.h:21
TGeoShape::kGeoPcon
@ kGeoPcon
Definition: TGeoShape.h:51
TGeoPcon::Divide
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this polycone shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition: TGeoPcon.cxx:627
TMath.h
int