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