Logo ROOT  
Reference Guide
TGeoArb8.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Andrei Gheata 31/01/02
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "TGeoArb8.h"
13
14#include <iostream>
15#include "TBuffer.h"
16#include "TGeoManager.h"
17#include "TGeoVolume.h"
18#include "TGeoMatrix.h"
19#include "TMath.h"
20
22
23/** \class TGeoArb8
24\ingroup Geometry_classes
25
26An arbitrary trapezoid with less than 8 vertices standing on
27two parallel planes perpendicular to Z axis. Parameters :
28 - dz - half length in Z;
29 - xy[8][2] - vector of (x,y) coordinates of vertices
30 - first four points (xy[i][j], i<4, j<2) are the (x,y)
31 coordinates of the vertices sitting on the -dz plane;
32 - last four points (xy[i][j], i>=4, j<2) are the (x,y)
33 coordinates of the vertices sitting on the +dz plane;
34
35The order of defining the vertices of an arb8 is the following :
36 - point 0 is connected with points 1,3,4
37 - point 1 is connected with points 0,2,5
38 - point 2 is connected with points 1,3,6
39 - point 3 is connected with points 0,2,7
40 - point 4 is connected with points 0,5,7
41 - point 5 is connected with points 1,4,6
42 - point 6 is connected with points 2,5,7
43 - point 7 is connected with points 3,4,6
44
45Points can be identical in order to create shapes with less than
468 vertices.
47
48Begin_Macro(source)
49{
50 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
51 new TGeoManager("arb8", "poza12");
52 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
53 TGeoMedium *med = new TGeoMedium("MED",1,mat);
54 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
55 gGeoManager->SetTopVolume(top);
56 TGeoArb8 *arb = new TGeoArb8(20);
57 arb->SetVertex(0,-30,-25);
58 arb->SetVertex(1,-25,25);
59 arb->SetVertex(2,5,25);
60 arb->SetVertex(3,25,-25);
61 arb->SetVertex(4,-28,-23);
62 arb->SetVertex(5,-23,27);
63 arb->SetVertex(6,-23,27);
64 arb->SetVertex(7,13,-27);
65 TGeoVolume *vol = new TGeoVolume("ARB8",arb,med);
66 vol->SetLineWidth(2);
67 top->AddNode(vol,1);
68 gGeoManager->CloseGeometry();
69 gGeoManager->SetNsegments(80);
70 top->Draw();
71 TView *view = gPad->GetView();
72 view->ShowAxis();
73}
74End_Macro
75*/
76
77/** \class TGeoGtra
78\ingroup Geometry_classes
79
80Gtra is a twisted trapezoid.
81i.e. one for which the faces perpendicular
82to z are trapezia and their centres are not the same x, y. It has 12
83parameters: the half length in z, the polar angles from the centre of
84the face at low z to that at high z, twist, H1 the half length in y at low z,
85LB1 the half length in x at low z and y low edge, LB2 the half length
86in x at low z and y high edge, TH1 the angle w.r.t. the y axis from the
87centre of low y edge to the centre of the high y edge, and H2, LB2,
88LH2, TH2, the corresponding quantities at high z.
89
90Begin_Macro(source)
91{
92 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
93 new TGeoManager("gtra", "poza11");
94 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
95 TGeoMedium *med = new TGeoMedium("MED",1,mat);
96 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
97 gGeoManager->SetTopVolume(top);
98 TGeoVolume *vol = gGeoManager->MakeGtra("Gtra",med, 30,15,30,30,20,10,15,0,20,10,15,0);
99 vol->SetLineWidth(2);
100 top->AddNode(vol,1);
101 gGeoManager->CloseGeometry();
102 gGeoManager->SetNsegments(80);
103 top->Draw();
104 TView *view = gPad->GetView();
105 view->ShowAxis();
106}
107End_Macro
108*/
109
110/** \class TGeoTrap
111\ingroup Geometry_classes
112
113TRAP is a general trapezoid, i.e. one for which the faces perpendicular
114to z are trapezia and their centres are not the same x, y. It has 11
115parameters: the half length in z, the polar angles from the centre of
116the face at low z to that at high z, H1 the half length in y at low z,
117LB1 the half length in x at low z and y low edge, LB2 the half length
118in x at low z and y high edge, TH1 the angle w.r.t. the y axis from the
119centre of low y edge to the centre of the high y edge, and H2, LB2,
120LH2, TH2, the corresponding quantities at high z.
121
122Begin_Macro(source)
123{
124 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
125 new TGeoManager("trap", "poza10");
126 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
127 TGeoMedium *med = new TGeoMedium("MED",1,mat);
128 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
129 gGeoManager->SetTopVolume(top);
130 TGeoVolume *vol = gGeoManager->MakeTrap("Trap",med, 30,15,30,20,10,15,0,20,10,15,0);
131 vol->SetLineWidth(2);
132 top->AddNode(vol,1);
133 gGeoManager->CloseGeometry();
134 gGeoManager->SetNsegments(80);
135 top->Draw();
136 TView *view = gPad->GetView();
137 view->ShowAxis();
138}
139End_Macro
140*/
141
142////////////////////////////////////////////////////////////////////////////////
143/// Default constructor.
144
146{
147 fDz = 0;
148 for (Int_t i=0; i<8; i++) {
149 fXY[i][0] = 0.0;
150 fXY[i][1] = 0.0;
151 }
153}
154
155////////////////////////////////////////////////////////////////////////////////
156/// Constructor. If the array of vertices is not null, this should be
157/// in the format : (x0, y0, x1, y1, ... , x7, y7)
158
160 :TGeoBBox(0,0,0)
161{
162 fDz = dz;
164 if (vertices) {
165 for (Int_t i=0; i<8; i++) {
166 fXY[i][0] = vertices[2*i];
167 fXY[i][1] = vertices[2*i+1];
168 }
169 ComputeTwist();
170 ComputeBBox();
171 } else {
172 for (Int_t i=0; i<8; i++) {
173 fXY[i][0] = 0.0;
174 fXY[i][1] = 0.0;
175 }
176 }
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// Named constructor. If the array of vertices is not null, this should be
181/// in the format : (x0, y0, x1, y1, ... , x7, y7)
182
183TGeoArb8::TGeoArb8(const char *name, Double_t dz, Double_t *vertices)
184 :TGeoBBox(name, 0,0,0)
185{
186 fDz = dz;
188 if (vertices) {
189 for (Int_t i=0; i<8; i++) {
190 fXY[i][0] = vertices[2*i];
191 fXY[i][1] = vertices[2*i+1];
192 }
193 ComputeTwist();
194 ComputeBBox();
195 } else {
196 for (Int_t i=0; i<8; i++) {
197 fXY[i][0] = 0.0;
198 fXY[i][1] = 0.0;
199 }
200 }
201}
202
203////////////////////////////////////////////////////////////////////////////////
204/// Destructor.
205
207{
208 if (fTwist) delete [] fTwist;
209}
210
211////////////////////////////////////////////////////////////////////////////////
212/// Copy twist values from source array
213
215{
216 if (twist) {
217 if (!fTwist) fTwist = new Double_t[4];
218 memcpy(fTwist, twist, 4*sizeof(Double_t));
219 } else if (fTwist) {
220 delete [] fTwist;
221 fTwist = nullptr;
222 }
223}
224
225////////////////////////////////////////////////////////////////////////////////
226/// Computes capacity of the shape in [length^3].
227
229{
230 Int_t i,j;
231 Double_t capacity = 0;
232 for (i=0; i<4; i++) {
233 j = (i+1)%4;
234 capacity += 0.25*fDz*((fXY[i][0]+fXY[i+4][0])*(fXY[j][1]+fXY[j+4][1]) -
235 (fXY[j][0]+fXY[j+4][0])*(fXY[i][1]+fXY[i+4][1]) +
236 (1./3)*((fXY[i+4][0]-fXY[i][0])*(fXY[j+4][1]-fXY[j][1]) -
237 (fXY[j][0]-fXY[j+4][0])*(fXY[i][1]-fXY[i+4][1])));
238 }
239 return TMath::Abs(capacity);
240}
241
242////////////////////////////////////////////////////////////////////////////////
243/// Computes bounding box for an Arb8 shape.
244
246{
248 xmin = xmax = fXY[0][0];
249 ymin = ymax = fXY[0][1];
250
251 for (Int_t i=1; i<8; i++) {
252 if (xmin>fXY[i][0]) xmin=fXY[i][0];
253 if (xmax<fXY[i][0]) xmax=fXY[i][0];
254 if (ymin>fXY[i][1]) ymin=fXY[i][1];
255 if (ymax<fXY[i][1]) ymax=fXY[i][1];
256 }
257 fDX = 0.5*(xmax-xmin);
258 fDY = 0.5*(ymax-ymin);
259 fDZ = fDz;
260 fOrigin[0] = 0.5*(xmax+xmin);
261 fOrigin[1] = 0.5*(ymax+ymin);
262 fOrigin[2] = 0;
264}
265
266////////////////////////////////////////////////////////////////////////////////
267/// Computes tangents of twist angles (angles between projections on XY plane
268/// of corresponding -dz +dz edges). Computes also if the vertices are defined
269/// clockwise or anti-clockwise.
270
272{
273 Double_t twist[4];
274 Bool_t twisted = kFALSE;
275 Double_t dx1, dy1, dx2, dy2;
276 Bool_t singleBottom = kTRUE;
277 Bool_t singleTop = kTRUE;
278 Int_t i;
279 for (i=0; i<4; i++) {
280 dx1 = fXY[(i+1)%4][0]-fXY[i][0];
281 dy1 = fXY[(i+1)%4][1]-fXY[i][1];
283 twist[i] = 0;
284 continue;
285 }
286 singleBottom = kFALSE;
287 dx2 = fXY[4+(i+1)%4][0]-fXY[4+i][0];
288 dy2 = fXY[4+(i+1)%4][1]-fXY[4+i][1];
290 twist[i] = 0;
291 continue;
292 }
293 singleTop = kFALSE;
294 twist[i] = dy1*dx2 - dx1*dy2;
295 if (TMath::Abs(twist[i])<TGeoShape::Tolerance()) {
296 twist[i] = 0;
297 continue;
298 }
299 twist[i] = TMath::Sign(1.,twist[i]);
300 twisted = kTRUE;
301 }
302
303 CopyTwist(twisted ? twist : nullptr);
304
305 if (singleBottom) {
306 for (i=0; i<4; i++) {
307 fXY[i][0] += 1.E-8*fXY[i+4][0];
308 fXY[i][1] += 1.E-8*fXY[i+4][1];
309 }
310 }
311 if (singleTop) {
312 for (i=0; i<4; i++) {
313 fXY[i+4][0] += 1.E-8*fXY[i][0];
314 fXY[i+4][1] += 1.E-8*fXY[i][1];
315 }
316 }
317 Double_t sum1 = 0.;
318 Double_t sum2 = 0.;
319 Int_t j;
320 for (i=0; i<4; i++) {
321 j = (i+1)%4;
322 sum1 += fXY[i][0]*fXY[j][1]-fXY[j][0]*fXY[i][1];
323 sum2 += fXY[i+4][0]*fXY[j+4][1]-fXY[j+4][0]*fXY[i+4][1];
324 }
325 if (sum1*sum2 < -TGeoShape::Tolerance()) {
326 Fatal("ComputeTwist", "Shape %s type Arb8: Lower/upper faces defined with opposite clockwise", GetName());
327 return;
328 }
329 if (sum1>TGeoShape::Tolerance()) {
330 Error("ComputeTwist", "Shape %s type Arb8: Vertices must be defined clockwise in XY planes. Re-ordering...", GetName());
331 Double_t xtemp, ytemp;
332 xtemp = fXY[1][0];
333 ytemp = fXY[1][1];
334 fXY[1][0] = fXY[3][0];
335 fXY[1][1] = fXY[3][1];
336 fXY[3][0] = xtemp;
337 fXY[3][1] = ytemp;
338 xtemp = fXY[5][0];
339 ytemp = fXY[5][1];
340 fXY[5][0] = fXY[7][0];
341 fXY[5][1] = fXY[7][1];
342 fXY[7][0] = xtemp;
343 fXY[7][1] = ytemp;
344 }
345 // Check for illegal crossings.
346 Bool_t illegal_cross = kFALSE;
347 illegal_cross = TGeoShape::IsSegCrossing(fXY[0][0],fXY[0][1],fXY[1][0],fXY[1][1],
348 fXY[2][0],fXY[2][1],fXY[3][0],fXY[3][1]);
349 if (!illegal_cross)
350 illegal_cross = TGeoShape::IsSegCrossing(fXY[4][0],fXY[4][1],fXY[5][0],fXY[5][1],
351 fXY[6][0],fXY[6][1],fXY[7][0],fXY[7][1]);
352 if (illegal_cross) {
353 Error("ComputeTwist", "Shape %s type Arb8: Malformed polygon with crossing opposite segments", GetName());
354 InspectShape();
355 }
356}
357
358////////////////////////////////////////////////////////////////////////////////
359/// Get twist for segment I in range [0,3]
360
362{
363 return (!fTwist || iseg<0 || iseg>3) ? 0. : fTwist[iseg];
364}
365
366////////////////////////////////////////////////////////////////////////////////
367/// Get index of the edge of the quadrilater represented by vert closest to point.
368/// If [P1,P2] is the closest segment and P is the point, the function returns the fraction of the
369/// projection of (P1P) over (P1P2). If projection of P is not in range [P1,P2] return -1.
370
371Double_t TGeoArb8::GetClosestEdge(const Double_t *point, Double_t *vert, Int_t &isegment) const
372{
373 isegment = 0;
374 Int_t isegmin = 0;
375 Int_t i1, i2;
376 Double_t p1[2], p2[2];
377 Double_t lsq, ssq, dx, dy, dpx, dpy, u;
378 Double_t umin = -1.;
379 Double_t safe=1E30;
380 for (i1=0; i1<4; i1++) {
382 isegment = isegmin;
383 return umin;
384 }
385 i2 = (i1+1)%4;
386 p1[0] = vert[2*i1];
387 p1[1] = vert[2*i1+1];
388 p2[0] = vert[2*i2];
389 p2[1] = vert[2*i2+1];
390 dx = p2[0] - p1[0];
391 dy = p2[1] - p1[1];
392 dpx = point[0] - p1[0];
393 dpy = point[1] - p1[1];
394 lsq = dx*dx + dy*dy;
395 // Current segment collapsed to a point
397 ssq = dpx*dpx + dpy*dpy;
398 if (ssq < safe) {
399 safe = ssq;
400 isegmin = i1;
401 umin = -1;
402 }
403 continue;
404 }
405 // Projection fraction
406 u = (dpx*dx + dpy*dy)/lsq;
407 // If projection of P is outside P1P2 limits, take the distance from P to the closest of P1 and P2
408 if (u>1) {
409 // Outside, closer to P2
410 dpx = point[0]-p2[0];
411 dpy = point[1]-p2[1];
412 u = -1.;
413 } else {
414 if (u>=0) {
415 // Projection inside
416 dpx -= u*dx;
417 dpy -= u*dy;
418 } else {
419 // Outside, closer to P1
420 u = -1.;
421 }
422 }
423 ssq = dpx*dpx + dpy*dpy;
424 if (ssq < safe) {
425 safe = ssq;
426 isegmin = i1;
427 umin = u;
428 }
429 }
430 isegment = isegmin;
431 // safe = TMath::Sqrt(safe);
432 return umin;
433}
434
435////////////////////////////////////////////////////////////////////////////////
436/// Compute normal to closest surface from POINT.
437
438void TGeoArb8::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
439{
440 Double_t safc;
441 Double_t x0, y0, z0, x1, y1, z1, x2, y2, z2;
442 Double_t ax, ay, az, bx, by, bz;
443 Double_t fn;
444 safc = fDz-TMath::Abs(point[2]);
445 if (safc<10.*TGeoShape::Tolerance()) {
446 memset(norm,0,3*sizeof(Double_t));
447 norm[2] = (dir[2]>0)?1:(-1);
448 return;
449 }
450 Double_t vert[8];
451 SetPlaneVertices(point[2], vert);
452 // Get the closest edge (point should be on this edge within tolerance)
453 Int_t iseg;
454 Double_t frac = GetClosestEdge(point, vert, iseg);
455 if (frac<0) frac = 0.;
456 Int_t jseg = (iseg+1)%4;
457 x0 = vert[2*iseg];
458 y0 = vert[2*iseg+1];
459 z0 = point[2];
460 x2 = vert[2*jseg];
461 y2 = vert[2*jseg+1];
462 z2 = point[2];
463 x0 += frac*(x2-x0);
464 y0 += frac*(y2-y0);
465 x1 = fXY[iseg+4][0];
466 y1 = fXY[iseg+4][1];
467 z1 = fDz;
468 x1 += frac*(fXY[jseg+4][0]-x1);
469 y1 += frac*(fXY[jseg+4][1]-y1);
470 ax = x1-x0;
471 ay = y1-y0;
472 az = z1-z0;
473 bx = x2-x0;
474 by = y2-y0;
475 bz = z2-z0;
476 // Cross product of the vector given by the section segment (that contains the point) at z=point[2]
477 // and the vector connecting the point projection to its correspondent on the top edge (hard to swallow, isn't it?)
478 norm[0] = ay*bz-az*by;
479 norm[1] = az*bx-ax*bz;
480 norm[2] = ax*by-ay*bx;
481 fn = TMath::Sqrt(norm[0]*norm[0]+norm[1]*norm[1]+norm[2]*norm[2]);
482 // If point is on one edge, fn may be very small and the normal does not make sense - avoid divzero
483 if (fn<1E-10) {
484 norm[0] = 1.;
485 norm[1] = 0.;
486 norm[2] = 0.;
487 } else {
488 norm[0] /= fn;
489 norm[1] /= fn;
490 norm[2] /= fn;
491 }
492 // Make sure the dot product of the normal and the direction is positive.
493 if (dir[0]>-2. && dir[0]*norm[0]+dir[1]*norm[1]+dir[2]*norm[2] < 0) {
494 norm[0] = -norm[0];
495 norm[1] = -norm[1];
496 norm[2] = -norm[2];
497 }
498}
499
500////////////////////////////////////////////////////////////////////////////////
501/// Test if point is inside this shape.
502
504{
505 // first check Z range
506 if (TMath::Abs(point[2]) > fDz) return kFALSE;
507 // compute intersection between Z plane containing point and the arb8
508 Double_t poly[8];
509 // memset(&poly[0], 0, 8*sizeof(Double_t));
510 // SetPlaneVertices(point[2], &poly[0]);
511 Double_t cf = 0.5*(fDz-point[2])/fDz;
512 Int_t i;
513 for (i=0; i<4; i++) {
514 poly[2*i] = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
515 poly[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
516 }
517 return InsidePolygon(point[0],point[1],poly);
518}
519
520////////////////////////////////////////////////////////////////////////////////
521/// Computes distance to plane ipl :
522/// - ipl=0 : points 0,4,1,5
523/// - ipl=1 : points 1,5,2,6
524/// - ipl=2 : points 2,6,3,7
525/// - ipl=3 : points 3,7,0,4
526
527Double_t TGeoArb8::DistToPlane(const Double_t *point, const Double_t *dir, Int_t ipl, Bool_t in) const
528{
529 Double_t xa,xb,xc,xd;
530 Double_t ya,yb,yc,yd;
531 Double_t eps = 10.*TGeoShape::Tolerance();
532 Double_t norm[3];
533 Double_t dirp[3];
534 Double_t ndotd = 0;
535 Int_t j = (ipl+1)%4;
536 xa=fXY[ipl][0];
537 ya=fXY[ipl][1];
538 xb=fXY[ipl+4][0];
539 yb=fXY[ipl+4][1];
540 xc=fXY[j][0];
541 yc=fXY[j][1];
542 xd=fXY[4+j][0];
543 yd=fXY[4+j][1];
544 Double_t dz2 =0.5/fDz;
545 Double_t tx1 =dz2*(xb-xa);
546 Double_t ty1 =dz2*(yb-ya);
547 Double_t tx2 =dz2*(xd-xc);
548 Double_t ty2 =dz2*(yd-yc);
549 Double_t dzp =fDz+point[2];
550 Double_t xs1 =xa+tx1*dzp;
551 Double_t ys1 =ya+ty1*dzp;
552 Double_t xs2 =xc+tx2*dzp;
553 Double_t ys2 =yc+ty2*dzp;
554 Double_t dxs =xs2-xs1;
555 Double_t dys =ys2-ys1;
556 Double_t dtx =tx2-tx1;
557 Double_t dty =ty2-ty1;
558 Double_t a=(dtx*dir[1]-dty*dir[0]+(tx1*ty2-tx2*ty1)*dir[2])*dir[2];
559 Double_t signa = TMath::Sign(1.,a);
560 Double_t b=dxs*dir[1]-dys*dir[0]+(dtx*point[1]-dty*point[0]+ty2*xs1-ty1*xs2
561 +tx1*ys2-tx2*ys1)*dir[2];
562 Double_t c=dxs*point[1]-dys*point[0]+xs1*ys2-xs2*ys1;
563 Double_t x1,x2,y1,y2,xp,yp,zi,s;
564 if (TMath::Abs(a)<eps) {
565 // Surface is planar
566 if (TMath::Abs(b)<eps) return TGeoShape::Big(); // Track parallel to surface
567 s=-c/b;
568 if (TMath::Abs(s)<1.E-6 && TMath::Abs(TMath::Abs(point[2])-fDz)>eps) {
569 memcpy(dirp,dir,3*sizeof(Double_t));
570 dirp[0] = -3;
571 // Compute normal pointing outside
572 ((TGeoArb8*)this)->ComputeNormal(point,dirp,norm);
573 ndotd = dir[0]*norm[0]+dir[1]*norm[1]+dir[2]*norm[2];
574 if (!in) ndotd*=-1.;
575 if (ndotd>0) {
576 s = TMath::Max(0.,s);
577 zi = (point[0]-xs1)*(point[0]-xs2)+(point[1]-ys1)*(point[1]-ys2);
578 if (zi<=0) return s;
579 }
580 return TGeoShape::Big();
581 }
582 if (s<0) return TGeoShape::Big();
583 } else {
584 Double_t d=b*b-4*a*c;
585 if (d<0) return TGeoShape::Big();
586 Double_t smin=0.5*(-b-signa*TMath::Sqrt(d))/a;
587 Double_t smax=0.5*(-b+signa*TMath::Sqrt(d))/a;
588 s = smin;
589 if (TMath::Abs(s)<1.E-6 && TMath::Abs(TMath::Abs(point[2])-fDz)>eps) {
590 memcpy(dirp,dir,3*sizeof(Double_t));
591 dirp[0] = -3;
592 // Compute normal pointing outside
593 ((TGeoArb8*)this)->ComputeNormal(point,dirp,norm);
594 ndotd = dir[0]*norm[0]+dir[1]*norm[1]+dir[2]*norm[2];
595 if (!in) ndotd*=-1.;
596 if (ndotd>0) return TMath::Max(0.,s);
597 s = 0.; // ignore
598 }
599 if (s>eps) {
600 // Check smin
601 zi=point[2]+s*dir[2];
602 if (TMath::Abs(zi)<fDz) {
603 x1=xs1+tx1*dir[2]*s;
604 x2=xs2+tx2*dir[2]*s;
605 xp=point[0]+s*dir[0];
606 y1=ys1+ty1*dir[2]*s;
607 y2=ys2+ty2*dir[2]*s;
608 yp=point[1]+s*dir[1];
609 zi = (xp-x1)*(xp-x2)+(yp-y1)*(yp-y2);
610 if (zi<=0) return s;
611 }
612 }
613 // Smin failed, try smax
614 s=smax;
615 if (TMath::Abs(s)<1.E-6 && TMath::Abs(TMath::Abs(point[2])-fDz)>eps) {
616 memcpy(dirp,dir,3*sizeof(Double_t));
617 dirp[0] = -3;
618 // Compute normal pointing outside
619 ((TGeoArb8*)this)->ComputeNormal(point,dirp,norm);
620 ndotd = dir[0]*norm[0]+dir[1]*norm[1]+dir[2]*norm[2];
621 if (!in) ndotd*=-1.;
622 if (ndotd>0) s = TMath::Max(0.,s);
623 else s = TGeoShape::Big();
624 return s;
625 }
626 }
627 if (s>eps) {
628 // Check smin
629 zi=point[2]+s*dir[2];
630 if (TMath::Abs(zi)<fDz) {
631 x1=xs1+tx1*dir[2]*s;
632 x2=xs2+tx2*dir[2]*s;
633 xp=point[0]+s*dir[0];
634 y1=ys1+ty1*dir[2]*s;
635 y2=ys2+ty2*dir[2]*s;
636 yp=point[1]+s*dir[1];
637 zi = (xp-x1)*(xp-x2)+(yp-y1)*(yp-y2);
638 if (zi<=0) return s;
639 }
640 }
641 return TGeoShape::Big();
642}
643
644////////////////////////////////////////////////////////////////////////////////
645/// Computes distance from outside point to surface of the shape.
646
647Double_t TGeoArb8::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t /*iact*/, Double_t step, Double_t * /*safe*/) const
648{
649 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
650 if (sdist>=step) return TGeoShape::Big();
652 // check Z planes
653 if (TMath::Abs(point[2])>fDz-1.E-8) {
654 Double_t pt[3];
655 if (point[2]*dir[2]<0) {
656 pt[2]=fDz*TMath::Sign(1.,point[2]);
657 snext = TMath::Max((pt[2]-point[2])/dir[2],0.);
658 for (Int_t j=0; j<2; j++) pt[j]=point[j]+snext*dir[j];
659 if (Contains(&pt[0])) return snext;
660 }
661 }
662 // check lateral faces
665 for (Int_t i=0; i<4; i++) {
666 dist = DistToPlane(point, dir, i, kFALSE);
667 if (dist<snext) snext = dist;
668 }
669 return snext;
670}
671
672////////////////////////////////////////////////////////////////////////////////
673/// Compute distance from inside point to surface of the shape.
674
675Double_t TGeoArb8::DistFromInside(const Double_t *point, const Double_t *dir, Int_t /*iact*/, Double_t /*step*/, Double_t * /*safe*/) const
676{
677 Int_t i;
678 Double_t distz = TGeoShape::Big();
679 Double_t distl = TGeoShape::Big();
681 Double_t pt[3] = {0.,0.,0.};
682 if (dir[2]<0) {
683 distz=(-fDz-point[2])/dir[2];
684 pt[2] = -fDz;
685 } else {
686 if (dir[2]>0) distz=(fDz-point[2])/dir[2];
687 pt[2] = fDz;
688 }
689 for (i=0; i<4; i++) {
690 dist=DistToPlane(point, dir, i, kTRUE);
691 if (dist<distl) distl = dist;
692 }
693 if (distz<distl) {
694 pt[0] = point[0]+distz*dir[0];
695 pt[1] = point[1]+distz*dir[1];
696 if (!Contains(pt)) distz = TGeoShape::Big();
697 }
698 dist = TMath::Min(distz, distl);
699 if (dist<0 || dist>1.E10) return 0.;
700 return dist;
701#ifdef OLDALGORITHM
702//#else
703// compute distance to plane ipl :
704// ipl=0 : points 0,4,1,5
705// ipl=1 : points 1,5,2,6
706// ipl=2 : points 2,6,3,7
707// ipl=3 : points 3,7,0,4
708 Double_t distmin;
709 Bool_t lateral_cross = kFALSE;
710 if (dir[2]<0) {
711 distmin=(-fDz-point[2])/dir[2];
712 } else {
713 if (dir[2]>0) distmin =(fDz-point[2])/dir[2];
714 else distmin = TGeoShape::Big();
715 }
716 Double_t dz2 =0.5/fDz;
717 Double_t xa,xb,xc,xd;
718 Double_t ya,yb,yc,yd;
719 Double_t eps = 100.*TGeoShape::Tolerance();
720 for (Int_t ipl=0;ipl<4;ipl++) {
721 Int_t j = (ipl+1)%4;
722 xa=fXY[ipl][0];
723 ya=fXY[ipl][1];
724 xb=fXY[ipl+4][0];
725 yb=fXY[ipl+4][1];
726 xc=fXY[j][0];
727 yc=fXY[j][1];
728 xd=fXY[4+j][0];
729 yd=fXY[4+j][1];
730
731 Double_t tx1 =dz2*(xb-xa);
732 Double_t ty1 =dz2*(yb-ya);
733 Double_t tx2 =dz2*(xd-xc);
734 Double_t ty2 =dz2*(yd-yc);
735 Double_t dzp =fDz+point[2];
736 Double_t xs1 =xa+tx1*dzp;
737 Double_t ys1 =ya+ty1*dzp;
738 Double_t xs2 =xc+tx2*dzp;
739 Double_t ys2 =yc+ty2*dzp;
740 Double_t dxs =xs2-xs1;
741 Double_t dys =ys2-ys1;
742 Double_t dtx =tx2-tx1;
743 Double_t dty =ty2-ty1;
744 Double_t a=(dtx*dir[1]-dty*dir[0]+(tx1*ty2-tx2*ty1)*dir[2])*dir[2];
745 Double_t b=dxs*dir[1]-dys*dir[0]+(dtx*point[1]-dty*point[0]+ty2*xs1-ty1*xs2
746 +tx1*ys2-tx2*ys1)*dir[2];
747 Double_t c=dxs*point[1]-dys*point[0]+xs1*ys2-xs2*ys1;
749 if (TMath::Abs(a)<eps) {
750 if (TMath::Abs(b)<eps) continue;
751 s=-c/b;
752 if (s>eps && s < distmin) {
753 distmin =s;
754 lateral_cross=kTRUE;
755 }
756 continue;
757 }
758 Double_t d=b*b-4*a*c;
759 if (d>=0.) {
760 if (a>0) s=0.5*(-b-TMath::Sqrt(d))/a;
761 else s=0.5*(-b+TMath::Sqrt(d))/a;
762 if (s>eps) {
763 if (s < distmin) {
764 distmin = s;
765 lateral_cross = kTRUE;
766 }
767 } else {
768 if (a>0) s=0.5*(-b+TMath::Sqrt(d))/a;
769 else s=0.5*(-b-TMath::Sqrt(d))/a;
770 if (s>eps && s < distmin) {
771 distmin =s;
772 lateral_cross = kTRUE;
773 }
774 }
775 }
776 }
777
778 if (!lateral_cross) {
779 // We have to make sure that track crosses the top or bottom.
780 if (distmin > 1.E10) return TGeoShape::Tolerance();
781 Double_t pt[2];
782 pt[0] = point[0]+distmin*dir[0];
783 pt[1] = point[1]+distmin*dir[1];
784 // Check if propagated point is in the polygon
785 Double_t poly[8];
786 Int_t i = 0;
787 if (dir[2]>0.) i=4;
788 for (Int_t j=0; j<4; j++) {
789 poly[2*j] = fXY[j+i][0];
790 poly[2*j+1] = fXY[j+i][1];
791 }
792 if (!InsidePolygon(pt[0],pt[1],poly)) return TGeoShape::Tolerance();
793 }
794 return distmin;
795#endif
796}
797
798////////////////////////////////////////////////////////////////////////////////
799/// Divide this shape along one axis.
800
801TGeoVolume *TGeoArb8::Divide(TGeoVolume *voldiv, const char * /*divname*/, Int_t /*iaxis*/, Int_t /*ndiv*/,
802 Double_t /*start*/, Double_t /*step*/)
803{
804 Error("Divide", "Division of an arbitrary trapezoid not implemented");
805 return voldiv;
806}
807
808////////////////////////////////////////////////////////////////////////////////
809/// Get shape range on a given axis.
810
812{
813 xlo = 0;
814 xhi = 0;
815 Double_t dx = 0;
816 if (iaxis==3) {
817 xlo = -fDz;
818 xhi = fDz;
819 dx = xhi-xlo;
820 return dx;
821 }
822 return dx;
823}
824
825////////////////////////////////////////////////////////////////////////////////
826/// Fill vector param[4] with the bounding cylinder parameters. The order
827/// is the following : Rmin, Rmax, Phi1, Phi2
828
830{
831 // first compute rmin/rmax
832 Double_t rmaxsq = 0;
833 Double_t rsq;
834 Int_t i;
835 for (i=0; i<8; i++) {
836 rsq = fXY[i][0]*fXY[i][0] + fXY[i][1]*fXY[i][1];
837 rmaxsq = TMath::Max(rsq, rmaxsq);
838 }
839 param[0] = 0.; // Rmin
840 param[1] = rmaxsq; // Rmax
841 param[2] = 0.; // Phi1
842 param[3] = 360.; // Phi2
843}
844
845////////////////////////////////////////////////////////////////////////////////
846/// Fills real parameters of a positioned box inside this arb8. Returns 0 if successful.
847
848Int_t TGeoArb8::GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const
849{
850 dx=dy=dz=0;
851 if (mat->IsRotation()) {
852 Error("GetFittingBox", "cannot handle parametrized rotated volumes");
853 return 1; // ### rotation not accepted ###
854 }
855 //--> translate the origin of the parametrized box to the frame of this box.
856 Double_t origin[3];
857 mat->LocalToMaster(parambox->GetOrigin(), origin);
858 if (!Contains(origin)) {
859 Error("GetFittingBox", "wrong matrix - parametrized box is outside this");
860 return 1; // ### wrong matrix ###
861 }
862 //--> now we have to get the valid range for all parametrized axis
863 Double_t dd[3];
864 dd[0] = parambox->GetDX();
865 dd[1] = parambox->GetDY();
866 dd[2] = parambox->GetDZ();
867 //-> check if Z range is fixed
868 if (dd[2]<0) {
869 dd[2] = TMath::Min(origin[2]+fDz, fDz-origin[2]);
870 if (dd[2]<0) {
871 Error("GetFittingBox", "wrong matrix");
872 return 1;
873 }
874 }
875 if (dd[0]>=0 && dd[1]>=0) {
876 dx = dd[0];
877 dy = dd[1];
878 dz = dd[2];
879 return 0;
880 }
881 //-> check now vertices at Z = origin[2] +/- dd[2]
882 Double_t upper[8];
883 Double_t lower[8];
884 SetPlaneVertices(origin[2]-dd[2], lower);
885 SetPlaneVertices(origin[2]+dd[2], upper);
886 for (Int_t iaxis=0; iaxis<2; iaxis++) {
887 if (dd[iaxis]>=0) continue;
888 Double_t ddmin = TGeoShape::Big();
889 for (Int_t ivert=0; ivert<4; ivert++) {
890 ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-lower[2*ivert+iaxis]));
891 ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-upper[2*ivert+iaxis]));
892 }
893 dd[iaxis] = ddmin;
894 }
895 dx = dd[0];
896 dy = dd[1];
897 dz = dd[2];
898 return 0;
899}
900
901////////////////////////////////////////////////////////////////////////////////
902/// Computes normal to plane defined by P1, P2 and P3
903
905{
906 Double_t cross = 0.;
907 Double_t v1[3], v2[3];
908 Int_t i;
909 for (i=0; i<3; i++) {
910 v1[i] = p2[i] - p1[i];
911 v2[i] = p3[i] - p1[i];
912 }
913 norm[0] = v1[1]*v2[2]-v1[2]*v2[1];
914 cross += norm[0]*norm[0];
915 norm[1] = v1[2]*v2[0]-v1[0]*v2[2];
916 cross += norm[1]*norm[1];
917 norm[2] = v1[0]*v2[1]-v1[1]*v2[0];
918 cross += norm[2]*norm[2];
919 if (TMath::Abs(cross) < TGeoShape::Tolerance()) return;
920 cross = 1./TMath::Sqrt(cross);
921 for (i=0; i<3; i++) norm[i] *= cross;
922}
923
924////////////////////////////////////////////////////////////////////////////////
925/// Fills array with n random points located on the surface of indexed facet.
926/// The output array must be provided with a length of minimum 3*npoints. Returns
927/// true if operation succeeded.
928/// Possible index values:
929/// - 0 - all facets together
930/// - 1 to 6 - facet index from bottom to top Z
931
932Bool_t TGeoArb8::GetPointsOnFacet(Int_t /*index*/, Int_t /*npoints*/, Double_t * /* array */) const
933{
934 return kFALSE;
935/*
936 if (index<0 || index>6) return kFALSE;
937 if (index==0) {
938 // Just generate same number of points on each facet
939 Int_t npts = npoints/6.;
940 Int_t count = 0;
941 for (Int_t ifacet=0; ifacet<6; ifacet++) {
942 if (GetPointsOnFacet(ifacet+1, npts, &array[3*count])) count += npts;
943 if (ifacet<5) npts = (npoints-count)/(5.-ifacet);
944 }
945 if (count>0) return kTRUE;
946 return kFALSE;
947 }
948 Double_t z, cf;
949 Double_t xmin=TGeoShape::Big();
950 Double_t xmax=-xmin;
951 Double_t ymin=TGeoShape::Big();
952 Double_t ymax=-ymin;
953 Double_t dy=0.;
954 Double_t poly[8];
955 Double_t point[2];
956 Int_t i;
957 if (index==1 || index==6) {
958 z = (index==1)?-fDz:fDz;
959 cf = 0.5*(fDz-z)/fDz;
960 for (i=0; i<4; i++) {
961 poly[2*i] = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
962 poly[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
963 xmin = TMath::Min(xmin, poly[2*i]);
964 xmax = TMath::Max(xmax, poly[2*i]);
965 ymin = TMath::Min(ymin, poly[2*i]);
966 ymax = TMath::Max(ymax, poly[2*i]);
967 }
968 }
969 Int_t nshoot = 0;
970 Int_t nmiss = 0;
971 for (i=0; i<npoints; i++) {
972 Double_t *point = &array[3*i];
973 switch (surfindex) {
974 case 1:
975 case 6:
976 while (nmiss<1000) {
977 point[0] = xmin + (xmax-xmin)*gRandom->Rndm();
978 point[1] = ymin + (ymax-ymin)*gRandom->Rndm();
979 }
980
981 return InsidePolygon(point[0],point[1],poly);
982*/
983}
984
985////////////////////////////////////////////////////////////////////////////////
986/// Finds if a point in XY plane is inside the polygon defines by PTS.
987
989{
990 Int_t i,j;
991 Double_t x1,y1,x2,y2;
992 Double_t cross;
993 for (i=0; i<4; i++) {
994 j = (i+1)%4;
995 x1 = pts[i<<1];
996 y1 = pts[(i<<1)+1];
997 x2 = pts[j<<1];
998 y2 = pts[(j<<1)+1];
999 cross = (x-x1)*(y2-y1)-(y-y1)*(x2-x1);
1000 if (cross<0) return kFALSE;
1001 }
1002 return kTRUE;
1003}
1004
1005////////////////////////////////////////////////////////////////////////////////
1006/// Prints shape parameters
1007
1009{
1010 printf("*** Shape %s: TGeoArb8 ***\n", GetName());
1011 if (IsTwisted()) printf(" = TWISTED\n");
1012 for (Int_t ip=0; ip<8; ip++) {
1013 printf(" point #%i : x=%11.5f y=%11.5f z=%11.5f\n",
1014 ip, fXY[ip][0], fXY[ip][1], fDz*((ip<4)?-1:1));
1015 }
1016 printf(" Bounding box:\n");
1018}
1019
1020////////////////////////////////////////////////////////////////////////////////
1021/// Computes the closest distance from given point to this shape.
1022
1024{
1025 Double_t safz = fDz-TMath::Abs(point[2]);
1026 if (!in) safz = -safz;
1027 Int_t iseg;
1028 Double_t safe = TGeoShape::Big();
1029 Double_t lsq, ssq, dx, dy, dpx, dpy, u;
1030 if (IsTwisted()) {
1031 if (!in) {
1032 if (!TGeoBBox::Contains(point)) return TGeoBBox::Safety(point,kFALSE);
1033 }
1034 // Point is also in the bounding box ;-(
1035 // Compute closest distance to any segment
1036 Double_t vert[8];
1037 Double_t *p1, *p2;
1038 Int_t isegmin=0;
1039 Double_t umin = 0.;
1040 SetPlaneVertices (point[2], vert);
1041 for (iseg=0; iseg<4; iseg++) {
1042 if (safe<TGeoShape::Tolerance()) return 0.;
1043 p1 = &vert[2*iseg];
1044 p2 = &vert[2*((iseg+1)%4)];
1045 dx = p2[0] - p1[0];
1046 dy = p2[1] - p1[1];
1047 dpx = point[0] - p1[0];
1048 dpy = point[1] - p1[1];
1049
1050 lsq = dx*dx + dy*dy;
1051 u = (dpx*dx + dpy*dy)/lsq;
1052 if (u>1) {
1053 dpx = point[0]-p2[0];
1054 dpy = point[1]-p2[1];
1055 } else {
1056 if (u>=0) {
1057 dpx -= u*dx;
1058 dpy -= u*dy;
1059 }
1060 }
1061 ssq = dpx*dpx + dpy*dpy;
1062 if (ssq < safe) {
1063 isegmin = iseg;
1064 umin = u;
1065 safe = ssq;
1066 }
1067 }
1068 if (umin<0) umin = 0.;
1069 if (umin>1) {
1070 isegmin = (isegmin+1)%4;
1071 umin = 0.;
1072 }
1073 Int_t i1 = isegmin;
1074 Int_t i2 = (isegmin+1)%4;
1075 Double_t dx1 = fXY[i2][0]-fXY[i1][0];
1076 Double_t dx2 = fXY[i2+4][0]-fXY[i1+4][0];
1077 Double_t dy1 = fXY[i2][1]-fXY[i1][1];
1078 Double_t dy2 = fXY[i2+4][1]-fXY[i1+4][1];
1079 dx = dx1 + umin*(dx2-dx1);
1080 dy = dy1 + umin*(dy2-dy1);
1081 safe *= 1.- 4.*fDz*fDz/(dx*dx+dy*dy+4.*fDz*fDz);
1082 safe = TMath::Sqrt(safe);
1083 if (in) return TMath::Min(safz,safe);
1084 return TMath::Max(safz,safe);
1085 }
1086
1087 Double_t saf[5];
1088 saf[0] = safz;
1089
1090 for (iseg=0; iseg<4; iseg++) saf[iseg+1] = SafetyToFace(point,iseg,in);
1091 if (in) safe = saf[TMath::LocMin(5, saf)];
1092 else safe = saf[TMath::LocMax(5, saf)];
1093 if (safe<0) return 0.;
1094 return safe;
1095}
1096
1097////////////////////////////////////////////////////////////////////////////////
1098/// Estimate safety to lateral plane defined by segment iseg in range [0,3]
1099/// Might be negative: plane seen only from inside.
1100
1102{
1103 Double_t vertices[12];
1104 Int_t ipln = (iseg+1)%4;
1105 // point 1
1106 vertices[0] = fXY[iseg][0];
1107 vertices[1] = fXY[iseg][1];
1108 vertices[2] = -fDz;
1109 // point 2
1110 vertices[3] = fXY[ipln][0];
1111 vertices[4] = fXY[ipln][1];
1112 vertices[5] = -fDz;
1113 // point 3
1114 vertices[6] = fXY[ipln+4][0];
1115 vertices[7] = fXY[ipln+4][1];
1116 vertices[8] = fDz;
1117 // point 4
1118 vertices[9] = fXY[iseg+4][0];
1119 vertices[10] = fXY[iseg+4][1];
1120 vertices[11] = fDz;
1121 Double_t safe;
1122 Double_t norm[3];
1123 Double_t *p1, *p2, *p3;
1124 p1 = &vertices[0];
1125 p2 = &vertices[9];
1126 p3 = &vertices[6];
1127 if (IsSamePoint(p2,p3)) {
1128 p3 = &vertices[3];
1129 if (IsSamePoint(p1,p3)) return -TGeoShape::Big(); // skip single segment
1130 }
1131 GetPlaneNormal(p1,p2,p3,norm);
1132 safe = (point[0]-p1[0])*norm[0]+(point[1]-p1[1])*norm[1]+(point[2]-p1[2])*norm[2];
1133 if (in) return (-safe);
1134 return safe;
1135}
1136
1137////////////////////////////////////////////////////////////////////////////////
1138/// Save a primitive as a C++ statement(s) on output stream "out".
1139
1140void TGeoArb8::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1141{
1143 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
1144 out << " dz = " << fDz << ";" << std::endl;
1145 out << " vert[0] = " << fXY[0][0] << ";" << std::endl;
1146 out << " vert[1] = " << fXY[0][1] << ";" << std::endl;
1147 out << " vert[2] = " << fXY[1][0] << ";" << std::endl;
1148 out << " vert[3] = " << fXY[1][1] << ";" << std::endl;
1149 out << " vert[4] = " << fXY[2][0] << ";" << std::endl;
1150 out << " vert[5] = " << fXY[2][1] << ";" << std::endl;
1151 out << " vert[6] = " << fXY[3][0] << ";" << std::endl;
1152 out << " vert[7] = " << fXY[3][1] << ";" << std::endl;
1153 out << " vert[8] = " << fXY[4][0] << ";" << std::endl;
1154 out << " vert[9] = " << fXY[4][1] << ";" << std::endl;
1155 out << " vert[10] = " << fXY[5][0] << ";" << std::endl;
1156 out << " vert[11] = " << fXY[5][1] << ";" << std::endl;
1157 out << " vert[12] = " << fXY[6][0] << ";" << std::endl;
1158 out << " vert[13] = " << fXY[6][1] << ";" << std::endl;
1159 out << " vert[14] = " << fXY[7][0] << ";" << std::endl;
1160 out << " vert[15] = " << fXY[7][1] << ";" << std::endl;
1161 out << " TGeoShape *" << GetPointerName() << " = new TGeoArb8(\"" << GetName() << "\", dz,vert);" << std::endl;
1163}
1164
1165////////////////////////////////////////////////////////////////////////////////
1166/// Computes intersection points between plane at zpl and non-horizontal edges.
1167
1169{
1170 Double_t cf = 0.5*(fDz-zpl)/fDz;
1171 for (Int_t i=0; i<4; i++) {
1172 vertices[2*i] = fXY[i+4][0]+cf*(fXY[i][0]-fXY[i+4][0]);
1173 vertices[2*i+1] = fXY[i+4][1]+cf*(fXY[i][1]-fXY[i+4][1]);
1174 }
1175}
1176
1177////////////////////////////////////////////////////////////////////////////////
1178/// Set all arb8 params in one step.
1179/// param[0] = dz
1180/// param[1] = x0
1181/// param[2] = y0
1182/// ...
1183
1185{
1186 fDz = param[0];
1187 for (Int_t i=0; i<8; i++) {
1188 fXY[i][0] = param[2*i+1];
1189 fXY[i][1] = param[2*i+2];
1190 }
1191 ComputeTwist();
1192 ComputeBBox();
1193}
1194
1195////////////////////////////////////////////////////////////////////////////////
1196/// Creates arb8 mesh points
1197
1199{
1200 for (Int_t i=0; i<8; i++) {
1201 points[3*i] = fXY[i][0];
1202 points[3*i+1] = fXY[i][1];
1203 points[3*i+2] = (i<4)?-fDz:fDz;
1204 }
1205}
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// Creates arb8 mesh points
1209
1211{
1212 for (Int_t i=0; i<8; i++) {
1213 points[3*i] = fXY[i][0];
1214 points[3*i+1] = fXY[i][1];
1215 points[3*i+2] = (i<4)?-fDz:fDz;
1216 }
1217}
1218
1219////////////////////////////////////////////////////////////////////////////////
1220/// Set values for a given vertex.
1221
1223{
1224 if (vnum<0 || vnum >7) {
1225 Error("SetVertex", "Invalid vertex number");
1226 return;
1227 }
1228 fXY[vnum][0] = x;
1229 fXY[vnum][1] = y;
1230 if (vnum == 7) {
1231 ComputeTwist();
1232 ComputeBBox();
1233 }
1234}
1235
1236////////////////////////////////////////////////////////////////////////////////
1237/// Fill size of this 3-D object
1238
1240{
1242}
1243
1244////////////////////////////////////////////////////////////////////////////////
1245/// Stream an object of class TGeoManager.
1246
1247void TGeoArb8::Streamer(TBuffer &R__b)
1248{
1249 if (R__b.IsReading()) {
1250 R__b.ReadClassBuffer(TGeoArb8::Class(), this);
1251 ComputeTwist();
1252 } else {
1253 R__b.WriteClassBuffer(TGeoArb8::Class(), this);
1254 }
1255}
1256
1257////////////////////////////////////////////////////////////////////////////////
1258/// Check the inside status for each of the points in the array.
1259/// Input: Array of point coordinates + vector size
1260/// Output: Array of Booleans for the inside of each point
1261
1262void TGeoArb8::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
1263{
1264 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
1265}
1266
1267////////////////////////////////////////////////////////////////////////////////
1268/// Compute the normal for an array o points so that norm.dot.dir is positive
1269/// Input: Arrays of point coordinates and directions + vector size
1270/// Output: Array of normal directions
1271
1272void TGeoArb8::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
1273{
1274 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
1275}
1276
1277////////////////////////////////////////////////////////////////////////////////
1278/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1279
1280void TGeoArb8::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1281{
1282 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
1283}
1284
1285////////////////////////////////////////////////////////////////////////////////
1286/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1287
1288void TGeoArb8::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1289{
1290 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
1291}
1292
1293////////////////////////////////////////////////////////////////////////////////
1294/// Compute safe distance from each of the points in the input array.
1295/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1296/// Output: Safety values
1297
1298void TGeoArb8::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
1299{
1300 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
1301}
1302
1304
1305////////////////////////////////////////////////////////////////////////////////
1306/// Default ctor
1307
1309{
1310 fDz = 0;
1311 fTheta = 0;
1312 fPhi = 0;
1313 fH1 = fH2 = fBl1 = fBl2 = fTl1 = fTl2 = fAlpha1 = fAlpha2 = 0;
1314}
1315
1316////////////////////////////////////////////////////////////////////////////////
1317/// Constructor providing just a range in Z, theta and phi.
1318
1320 :TGeoArb8("", 0, 0)
1321{
1322 fDz = dz;
1323 fTheta = theta;
1324 fPhi = phi;
1325 fH1 = fH2 = fBl1 = fBl2 = fTl1 = fTl2 = fAlpha1 = fAlpha2 = 0;
1326}
1327
1328////////////////////////////////////////////////////////////////////////////////
1329/// Normal constructor.
1330
1332 Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2,
1333 Double_t tl2, Double_t alpha2)
1334 :TGeoArb8("", 0, 0)
1335{
1336 fDz = dz;
1337 fTheta = theta;
1338 fPhi = phi;
1339 fH1 = h1;
1340 fH2 = h2;
1341 fBl1 = bl1;
1342 fBl2 = bl2;
1343 fTl1 = tl1;
1344 fTl2 = tl2;
1345 fAlpha1 = alpha1;
1346 fAlpha2 = alpha2;
1349 Double_t ta1 = TMath::Tan(alpha1*TMath::DegToRad());
1350 Double_t ta2 = TMath::Tan(alpha2*TMath::DegToRad());
1351 fXY[0][0] = -dz*tx-h1*ta1-bl1; fXY[0][1] = -dz*ty-h1;
1352 fXY[1][0] = -dz*tx+h1*ta1-tl1; fXY[1][1] = -dz*ty+h1;
1353 fXY[2][0] = -dz*tx+h1*ta1+tl1; fXY[2][1] = -dz*ty+h1;
1354 fXY[3][0] = -dz*tx-h1*ta1+bl1; fXY[3][1] = -dz*ty-h1;
1355 fXY[4][0] = dz*tx-h2*ta2-bl2; fXY[4][1] = dz*ty-h2;
1356 fXY[5][0] = dz*tx+h2*ta2-tl2; fXY[5][1] = dz*ty+h2;
1357 fXY[6][0] = dz*tx+h2*ta2+tl2; fXY[6][1] = dz*ty+h2;
1358 fXY[7][0] = dz*tx-h2*ta2+bl2; fXY[7][1] = dz*ty-h2;
1359 ComputeTwist();
1360 if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) ||
1361 (h2<0) || (bl2<0) || (tl2<0)) {
1363 }
1364 else TGeoArb8::ComputeBBox();
1365}
1366
1367////////////////////////////////////////////////////////////////////////////////
1368/// Constructor with name.
1369
1371 Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2,
1372 Double_t tl2, Double_t alpha2)
1373 :TGeoArb8(name, 0, 0)
1374{
1375 SetName(name);
1376 fDz = dz;
1377 fTheta = theta;
1378 fPhi = phi;
1379 fH1 = h1;
1380 fH2 = h2;
1381 fBl1 = bl1;
1382 fBl2 = bl2;
1383 fTl1 = tl1;
1384 fTl2 = tl2;
1385 fAlpha1 = alpha1;
1386 fAlpha2 = alpha2;
1387 for (Int_t i=0; i<8; i++) {
1388 fXY[i][0] = 0.0;
1389 fXY[i][1] = 0.0;
1390 }
1393 Double_t ta1 = TMath::Tan(alpha1*TMath::DegToRad());
1394 Double_t ta2 = TMath::Tan(alpha2*TMath::DegToRad());
1395 fXY[0][0] = -dz*tx-h1*ta1-bl1; fXY[0][1] = -dz*ty-h1;
1396 fXY[1][0] = -dz*tx+h1*ta1-tl1; fXY[1][1] = -dz*ty+h1;
1397 fXY[2][0] = -dz*tx+h1*ta1+tl1; fXY[2][1] = -dz*ty+h1;
1398 fXY[3][0] = -dz*tx-h1*ta1+bl1; fXY[3][1] = -dz*ty-h1;
1399 fXY[4][0] = dz*tx-h2*ta2-bl2; fXY[4][1] = dz*ty-h2;
1400 fXY[5][0] = dz*tx+h2*ta2-tl2; fXY[5][1] = dz*ty+h2;
1401 fXY[6][0] = dz*tx+h2*ta2+tl2; fXY[6][1] = dz*ty+h2;
1402 fXY[7][0] = dz*tx-h2*ta2+bl2; fXY[7][1] = dz*ty-h2;
1403 ComputeTwist();
1404 if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) ||
1405 (h2<0) || (bl2<0) || (tl2<0)) {
1407 }
1408 else TGeoArb8::ComputeBBox();
1409}
1410
1411////////////////////////////////////////////////////////////////////////////////
1412/// Destructor.
1413
1415{
1416}
1417
1418////////////////////////////////////////////////////////////////////////////////
1419/// Compute distance from inside point to surface of the trapezoid
1420
1421Double_t TGeoTrap::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1422{
1423 if (iact<3 && safe) {
1424 // compute safe distance
1425 *safe = Safety(point, kTRUE);
1426 if (iact==0) return TGeoShape::Big();
1427 if (iact==1 && step<*safe) return TGeoShape::Big();
1428 }
1429 // compute distance to get outside this shape
1430 // return TGeoArb8::DistFromInside(point, dir, iact, step, safe);
1431
1432// compute distance to plane ipl :
1433// ipl=0 : points 0,4,1,5
1434// ipl=1 : points 1,5,2,6
1435// ipl=2 : points 2,6,3,7
1436// ipl=3 : points 3,7,0,4
1437 Double_t distmin;
1438 if (dir[2]<0) {
1439 distmin=(-fDz-point[2])/dir[2];
1440 } else {
1441 if (dir[2]>0) distmin =(fDz-point[2])/dir[2];
1442 else distmin = TGeoShape::Big();
1443 }
1444 Double_t xa,xb,xc;
1445 Double_t ya,yb,yc;
1446 for (Int_t ipl=0;ipl<4;ipl++) {
1447 Int_t j = (ipl+1)%4;
1448 xa=fXY[ipl][0];
1449 ya=fXY[ipl][1];
1450 xb=fXY[ipl+4][0];
1451 yb=fXY[ipl+4][1];
1452 xc=fXY[j][0];
1453 yc=fXY[j][1];
1454 Double_t ax,ay,az;
1455 ax = xb-xa;
1456 ay = yb-ya;
1457 az = 2.*fDz;
1458 Double_t bx,by;
1459 bx = xc-xa;
1460 by = yc-ya;
1461 Double_t ddotn = -dir[0]*az*by + dir[1]*az*bx+dir[2]*(ax*by-ay*bx);
1462 if (ddotn<=0) continue; // entering
1463 Double_t saf = -(point[0]-xa)*az*by + (point[1]-ya)*az*bx + (point[2]+fDz)*(ax*by-ay*bx);
1464 if (saf>=0.0) return 0.0;
1465 Double_t s = -saf/ddotn;
1466 if (s<distmin) distmin=s;
1467 }
1468 return distmin;
1469}
1470
1471////////////////////////////////////////////////////////////////////////////////
1472/// Compute distance from outside point to surface of the trapezoid
1473
1474Double_t TGeoTrap::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1475{
1476 if (iact<3 && safe) {
1477 // compute safe distance
1478 *safe = Safety(point, kFALSE);
1479 if (iact==0) return TGeoShape::Big();
1480 if (iact==1 && step<*safe) return TGeoShape::Big();
1481 }
1482 // Check if the bounding box is crossed within the requested distance
1483 Double_t sdist = TGeoBBox::DistFromOutside(point,dir, fDX, fDY, fDZ, fOrigin, step);
1484 if (sdist>=step) return TGeoShape::Big();
1485 // compute distance to get outside this shape
1486 Bool_t in = kTRUE;
1487 Double_t pts[8];
1488 Double_t xnew, ynew, znew;
1489 Int_t i,j;
1490 if (point[2]<-fDz+TGeoShape::Tolerance()) {
1491 if (dir[2] < TGeoShape::Tolerance()) return TGeoShape::Big();
1492 in = kFALSE;
1493 Double_t snxt = -(fDz+point[2])/dir[2];
1494 xnew = point[0] + snxt*dir[0];
1495 ynew = point[1] + snxt*dir[1];
1496 for (i=0;i<4;i++) {
1497 j = i<<1;
1498 pts[j] = fXY[i][0];
1499 pts[j+1] = fXY[i][1];
1500 }
1501 if (InsidePolygon(xnew,ynew,pts)) return snxt;
1502 } else if (point[2]>fDz-TGeoShape::Tolerance()) {
1503 if (dir[2] > -TGeoShape::Tolerance()) return TGeoShape::Big();
1504 in = kFALSE;
1505 Double_t snxt = (fDz-point[2])/dir[2];
1506 xnew = point[0] + snxt*dir[0];
1507 ynew = point[1] + snxt*dir[1];
1508 for (i=0;i<4;i++) {
1509 j = i<<1;
1510 pts[j] = fXY[i+4][0];
1511 pts[j+1] = fXY[i+4][1];
1512 }
1513 if (InsidePolygon(xnew,ynew,pts)) return snxt;
1514 }
1515 // check lateral faces
1516 Double_t dz2 =0.5/fDz;
1517 Double_t xa,xb,xc,xd;
1518 Double_t ya,yb,yc,yd;
1519 Double_t ax,ay,az;
1520 Double_t bx,by;
1521 Double_t ddotn, saf;
1522 Double_t safmin = TGeoShape::Big();
1523 Bool_t exiting = kFALSE;
1524 for (i=0; i<4; i++) {
1525 j = (i+1)%4;
1526 xa=fXY[i][0];
1527 ya=fXY[i][1];
1528 xb=fXY[i+4][0];
1529 yb=fXY[i+4][1];
1530 xc=fXY[j][0];
1531 yc=fXY[j][1];
1532 xd=fXY[4+j][0];
1533 yd=fXY[4+j][1];
1534 ax = xb-xa;
1535 ay = yb-ya;
1536 az = 2.*fDz;
1537 bx = xc-xa;
1538 by = yc-ya;
1539 ddotn = -dir[0]*az*by + dir[1]*az*bx+dir[2]*(ax*by-ay*bx);
1540 saf = (point[0]-xa)*az*by - (point[1]-ya)*az*bx - (point[2]+fDz)*(ax*by-ay*bx);
1541
1542 if (saf<=0) {
1543 // face visible from point outside
1544 in = kFALSE;
1545 if (ddotn>=0) return TGeoShape::Big();
1546 Double_t snxt = saf/ddotn;
1547 znew = point[2]+snxt*dir[2];
1548 if (TMath::Abs(znew)<=fDz) {
1549 xnew = point[0]+snxt*dir[0];
1550 ynew = point[1]+snxt*dir[1];
1551 Double_t tx1 =dz2*(xb-xa);
1552 Double_t ty1 =dz2*(yb-ya);
1553 Double_t tx2 =dz2*(xd-xc);
1554 Double_t ty2 =dz2*(yd-yc);
1555 Double_t dzp =fDz+znew;
1556 Double_t xs1 =xa+tx1*dzp;
1557 Double_t ys1 =ya+ty1*dzp;
1558 Double_t xs2 =xc+tx2*dzp;
1559 Double_t ys2 =yc+ty2*dzp;
1560 if (TMath::Abs(xs1-xs2)>TMath::Abs(ys1-ys2)) {
1561 if ((xnew-xs1)*(xs2-xnew)>=0) return snxt;
1562 } else {
1563 if ((ynew-ys1)*(ys2-ynew)>=0) return snxt;
1564 }
1565 }
1566 } else {
1567 if (saf<safmin) {
1568 safmin = saf;
1569 if (ddotn>=0) exiting = kTRUE;
1570 else exiting = kFALSE;
1571 }
1572 }
1573 }
1574 // Check also Z boundaries (point may be inside and close to Z) - Christian Hammann
1575 saf = fDz - TMath::Abs(point[2]);
1576 if (saf>0 && saf<safmin) exiting = (point[2]*dir[2] > 0)?kTRUE:kFALSE;
1577 if (!in) return TGeoShape::Big();
1578 if (exiting) return TGeoShape::Big();
1579 return 0.0;
1580}
1581
1582////////////////////////////////////////////////////////////////////////////////
1583/// Divide this trapezoid shape belonging to volume "voldiv" into ndiv volumes
1584/// called divname, from start position with the given step. Only Z divisions
1585/// are supported. For Z divisions just return the pointer to the volume to be
1586/// divided. In case a wrong division axis is supplied, returns pointer to
1587/// volume that was divided.
1588
1589TGeoVolume *TGeoTrap::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
1590 Double_t start, Double_t step)
1591{
1592 TGeoShape *shape; //--- shape to be created
1593 TGeoVolume *vol; //--- division volume to be created
1594 TGeoVolumeMulti *vmulti; //--- generic divided volume
1595 TGeoPatternFinder *finder; //--- finder to be attached
1596 TString opt = ""; //--- option to be attached
1597 if (iaxis!=3) {
1598 Error("Divide", "cannot divide trapezoids on other axis than Z");
1599 return 0;
1600 }
1601 Double_t end = start+ndiv*step;
1602 Double_t points_lo[8];
1603 Double_t points_hi[8];
1604 finder = new TGeoPatternTrapZ(voldiv, ndiv, start, end);
1605 voldiv->SetFinder(finder);
1606 finder->SetDivIndex(voldiv->GetNdaughters());
1607 opt = "Z";
1608 vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
1609 Double_t txz = ((TGeoPatternTrapZ*)finder)->GetTxz();
1610 Double_t tyz = ((TGeoPatternTrapZ*)finder)->GetTyz();
1611 Double_t zmin, zmax, ox,oy,oz;
1612 for (Int_t idiv=0; idiv<ndiv; idiv++) {
1613 zmin = start+idiv*step;
1614 zmax = start+(idiv+1)*step;
1615 oz = start+idiv*step+step/2;
1616 ox = oz*txz;
1617 oy = oz*tyz;
1618 SetPlaneVertices(zmin, &points_lo[0]);
1619 SetPlaneVertices(zmax, &points_hi[0]);
1620 shape = new TGeoTrap(step/2, fTheta, fPhi);
1621 for (Int_t vert1=0; vert1<4; vert1++)
1622 ((TGeoArb8*)shape)->SetVertex(vert1, points_lo[2*vert1]-ox, points_lo[2*vert1+1]-oy);
1623 for (Int_t vert2=0; vert2<4; vert2++)
1624 ((TGeoArb8*)shape)->SetVertex(vert2+4, points_hi[2*vert2]-ox, points_hi[2*vert2+1]-oy);
1625 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
1626 vmulti->AddVolume(vol);
1627 voldiv->AddNodeOffset(vol, idiv, oz, opt.Data());
1628 ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
1629 }
1630 return vmulti;
1631}
1632
1633////////////////////////////////////////////////////////////////////////////////
1634/// In case shape has some negative parameters, these have to be computed
1635/// in order to fit the mother.
1636
1638{
1639 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
1640 if (mother->IsRunTimeShape()) {
1641 Error("GetMakeRuntimeShape", "invalid mother");
1642 return 0;
1643 }
1644 Double_t dz, h1, bl1, tl1, h2, bl2, tl2;
1645 if (fDz<0) dz=((TGeoTrap*)mother)->GetDz();
1646 else dz=fDz;
1647
1648 if (fH1<0) h1 = ((TGeoTrap*)mother)->GetH1();
1649 else h1 = fH1;
1650
1651 if (fH2<0) h2 = ((TGeoTrap*)mother)->GetH2();
1652 else h2 = fH2;
1653
1654 if (fBl1<0) bl1 = ((TGeoTrap*)mother)->GetBl1();
1655 else bl1 = fBl1;
1656
1657 if (fBl2<0) bl2 = ((TGeoTrap*)mother)->GetBl2();
1658 else bl2 = fBl2;
1659
1660 if (fTl1<0) tl1 = ((TGeoTrap*)mother)->GetTl1();
1661 else tl1 = fTl1;
1662
1663 if (fTl2<0) tl2 = ((TGeoTrap*)mother)->GetTl2();
1664 else tl2 = fTl2;
1665
1666 return (new TGeoTrap(dz, fTheta, fPhi, h1, bl1, tl1, fAlpha1, h2, bl2, tl2, fAlpha2));
1667}
1668
1669////////////////////////////////////////////////////////////////////////////////
1670/// Computes the closest distance from given point to this shape.
1671
1673{
1674 Double_t saf[5];
1675 Double_t norm[3]; // normal to current facette
1676 Int_t i, j; // current facette index
1677 Double_t x0, y0, z0=-fDz, x1, y1, z1=fDz, x2, y2;
1678 Double_t ax, ay, az=z1-z0, bx, by;
1679 Double_t fn, safe;
1680 //---> compute safety for lateral planes
1681 for (i=0; i<4; i++) {
1682 if (in) saf[i] = TGeoShape::Big();
1683 else saf[i] = 0.;
1684 x0 = fXY[i][0];
1685 y0 = fXY[i][1];
1686 x1 = fXY[i+4][0];
1687 y1 = fXY[i+4][1];
1688 ax = x1-x0;
1689 ay = y1-y0;
1690 j = (i+1)%4;
1691 x2 = fXY[j][0];
1692 y2 = fXY[j][1];
1693 bx = x2-x0;
1694 by = y2-y0;
1696 x2 = fXY[4+j][0];
1697 y2 = fXY[4+j][1];
1698 bx = x2-x1;
1699 by = y2-y1;
1701 }
1702 norm[0] = -az*by;
1703 norm[1] = az*bx;
1704 norm[2] = ax*by-ay*bx;
1705 fn = TMath::Sqrt(norm[0]*norm[0]+norm[1]*norm[1]+norm[2]*norm[2]);
1706 if (fn<1E-10) continue;
1707 saf[i] = (x0-point[0])*norm[0]+(y0-point[1])*norm[1]+(-fDz-point[2])*norm[2];
1708 if (in) {
1709 saf[i]=TMath::Abs(saf[i])/fn; // they should be all positive anyway
1710 } else {
1711 saf[i] = -saf[i]/fn; // only negative values are interesting
1712 }
1713 }
1714 saf[4] = fDz-TMath::Abs(point[2]);
1715 if (in) {
1716 safe = saf[0];
1717 for (j=1;j<5;j++)
1718 if (saf[j] < safe)
1719 safe = saf[j];
1720 } else {
1721 saf[4]=-saf[4];
1722 safe = saf[0];
1723 for (j=1;j<5;j++)
1724 if (saf[j] > safe)
1725 safe = saf[j];
1726 }
1727 return safe;
1728}
1729
1730////////////////////////////////////////////////////////////////////////////////
1731/// Save a primitive as a C++ statement(s) on output stream "out".
1732
1733void TGeoTrap::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1734{
1736 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
1737 out << " dz = " << fDz << ";" << std::endl;
1738 out << " theta = " << fTheta << ";" << std::endl;
1739 out << " phi = " << fPhi << ";" << std::endl;
1740 out << " h1 = " << fH1<< ";" << std::endl;
1741 out << " bl1 = " << fBl1<< ";" << std::endl;
1742 out << " tl1 = " << fTl1<< ";" << std::endl;
1743 out << " alpha1 = " << fAlpha1 << ";" << std::endl;
1744 out << " h2 = " << fH2 << ";" << std::endl;
1745 out << " bl2 = " << fBl2<< ";" << std::endl;
1746 out << " tl2 = " << fTl2<< ";" << std::endl;
1747 out << " alpha2 = " << fAlpha2 << ";" << std::endl;
1748 out << " TGeoShape *" << GetPointerName() << " = new TGeoTrap(\"" << GetName() << "\", dz,theta,phi,h1,bl1,tl1,alpha1,h2,bl2,tl2,alpha2);" << std::endl;
1750}
1751
1752////////////////////////////////////////////////////////////////////////////////
1753/// Set all arb8 params in one step.
1754/// - param[0] = dz
1755/// - param[1] = theta
1756/// - param[2] = phi
1757/// - param[3] = h1
1758/// - param[4] = bl1
1759/// - param[5] = tl1
1760/// - param[6] = alpha1
1761/// - param[7] = h2
1762/// - param[8] = bl2
1763/// - param[9] = tl2
1764/// - param[10] = alpha2
1765
1767{
1768 fDz = param[0];
1769 fTheta = param[1];
1770 fPhi = param[2];
1771 fH1 = param[3];
1772 fH2 = param[7];
1773 fBl1 = param[4];
1774 fBl2 = param[8];
1775 fTl1 = param[5];
1776 fTl2 = param[9];
1777 fAlpha1 = param[6];
1778 fAlpha2 = param[10];
1783 fXY[0][0] = -fDz*tx-fH1*ta1-fBl1; fXY[0][1] = -fDz*ty-fH1;
1784 fXY[1][0] = -fDz*tx+fH1*ta1-fTl1; fXY[1][1] = -fDz*ty+fH1;
1785 fXY[2][0] = -fDz*tx+fH1*ta1+fTl1; fXY[2][1] = -fDz*ty+fH1;
1786 fXY[3][0] = -fDz*tx-fH1*ta1+fBl1; fXY[3][1] = -fDz*ty-fH1;
1787 fXY[4][0] = fDz*tx-fH2*ta2-fBl2; fXY[4][1] = fDz*ty-fH2;
1788 fXY[5][0] = fDz*tx+fH2*ta2-fTl2; fXY[5][1] = fDz*ty+fH2;
1789 fXY[6][0] = fDz*tx+fH2*ta2+fTl2; fXY[6][1] = fDz*ty+fH2;
1790 fXY[7][0] = fDz*tx-fH2*ta2+fBl2; fXY[7][1] = fDz*ty-fH2;
1791 ComputeTwist();
1792 if ((fDz<0) || (fH1<0) || (fBl1<0) || (fTl1<0) ||
1793 (fH2<0) || (fBl2<0) || (fTl2<0)) {
1795 }
1796 else TGeoArb8::ComputeBBox();
1797}
1798
1799////////////////////////////////////////////////////////////////////////////////
1800/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1801
1802void TGeoTrap::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1803{
1804 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
1805}
1806
1807////////////////////////////////////////////////////////////////////////////////
1808/// Compute distance from array of input points having directions specified by dirs. Store output in dists
1809
1810void TGeoTrap::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
1811{
1812 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
1813}
1814
1815////////////////////////////////////////////////////////////////////////////////
1816/// Compute safe distance from each of the points in the input array.
1817/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
1818/// Output: Safety values
1819
1820void TGeoTrap::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
1821{
1822 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
1823}
1824
1826
1827////////////////////////////////////////////////////////////////////////////////
1828/// Default ctor
1829
1831{
1832 fTwistAngle = 0;
1833
1834}
1835
1836////////////////////////////////////////////////////////////////////////////////
1837/// Constructor.
1838
1840 Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2,
1841 Double_t tl2, Double_t alpha2)
1842 :TGeoTrap(dz, theta, phi, h1, bl1, tl1, alpha1, h2, bl2, tl2, alpha2)
1843{
1844 fTwistAngle = twist;
1845 Double_t x,y;
1846 Double_t th = theta*TMath::DegToRad();
1847 Double_t ph = phi*TMath::DegToRad();
1848 // Coordinates of the center of the bottom face
1849 Double_t xc = -dz*TMath::Sin(th)*TMath::Cos(ph);
1850 Double_t yc = -dz*TMath::Sin(th)*TMath::Sin(ph);
1851
1852 Int_t i;
1853
1854 for (i=0; i<4; i++) {
1855 x = fXY[i][0] - xc;
1856 y = fXY[i][1] - yc;
1857 fXY[i][0] = x*TMath::Cos(-0.5*twist*TMath::DegToRad()) + y*TMath::Sin(-0.5*twist*TMath::DegToRad()) + xc;
1858 fXY[i][1] = -x*TMath::Sin(-0.5*twist*TMath::DegToRad()) + y*TMath::Cos(-0.5*twist*TMath::DegToRad()) + yc;
1859 }
1860 // Coordinates of the center of the top face
1861 xc = -xc;
1862 yc = -yc;
1863 for (i=4; i<8; i++) {
1864 x = fXY[i][0] - xc;
1865 y = fXY[i][1] - yc;
1866 fXY[i][0] = x*TMath::Cos(0.5*twist*TMath::DegToRad()) + y*TMath::Sin(0.5*twist*TMath::DegToRad()) + xc;
1867 fXY[i][1] = -x*TMath::Sin(0.5*twist*TMath::DegToRad()) + y*TMath::Cos(0.5*twist*TMath::DegToRad()) + yc;
1868 }
1869 ComputeTwist();
1870 if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) ||
1871 (h2<0) || (bl2<0) || (tl2<0)) SetShapeBit(kGeoRunTimeShape);
1872 else TGeoArb8::ComputeBBox();
1873}
1874
1875////////////////////////////////////////////////////////////////////////////////
1876/// Constructor providing the name of the shape.
1877
1879 Double_t bl1, Double_t tl1, Double_t alpha1, Double_t h2, Double_t bl2,
1880 Double_t tl2, Double_t alpha2)
1881 :TGeoTrap(name, dz, theta, phi, h1, bl1, tl1, alpha1, h2, bl2, tl2, alpha2)
1882{
1883 fTwistAngle = twist;
1884 Double_t x,y;
1885 Double_t th = theta*TMath::DegToRad();
1886 Double_t ph = phi*TMath::DegToRad();
1887 // Coordinates of the center of the bottom face
1888 Double_t xc = -dz*TMath::Sin(th)*TMath::Cos(ph);
1889 Double_t yc = -dz*TMath::Sin(th)*TMath::Sin(ph);
1890
1891 Int_t i;
1892
1893 for (i=0; i<4; i++) {
1894 x = fXY[i][0] - xc;
1895 y = fXY[i][1] - yc;
1896 fXY[i][0] = x*TMath::Cos(-0.5*twist*TMath::DegToRad()) + y*TMath::Sin(-0.5*twist*TMath::DegToRad()) + xc;
1897 fXY[i][1] = -x*TMath::Sin(-0.5*twist*TMath::DegToRad()) + y*TMath::Cos(-0.5*twist*TMath::DegToRad()) + yc;
1898 }
1899 // Coordinates of the center of the top face
1900 xc = -xc;
1901 yc = -yc;
1902 for (i=4; i<8; i++) {
1903 x = fXY[i][0] - xc;
1904 y = fXY[i][1] - yc;
1905 fXY[i][0] = x*TMath::Cos(0.5*twist*TMath::DegToRad()) + y*TMath::Sin(0.5*twist*TMath::DegToRad()) + xc;
1906 fXY[i][1] = -x*TMath::Sin(0.5*twist*TMath::DegToRad()) + y*TMath::Cos(0.5*twist*TMath::DegToRad()) + yc;
1907 }
1908 ComputeTwist();
1909 if ((dz<0) || (h1<0) || (bl1<0) || (tl1<0) ||
1910 (h2<0) || (bl2<0) || (tl2<0)) SetShapeBit(kGeoRunTimeShape);
1911 else TGeoArb8::ComputeBBox();
1912}
1913
1914////////////////////////////////////////////////////////////////////////////////
1915/// Destructor.
1916
1918{
1919}
1920
1921////////////////////////////////////////////////////////////////////////////////
1922/// Compute distance from inside point to surface of the shape.
1923
1924Double_t TGeoGtra::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1925{
1926 if (iact<3 && safe) {
1927 // compute safe distance
1928 *safe = Safety(point, kTRUE);
1929 if (iact==0) return TGeoShape::Big();
1930 if (iact==1 && step<*safe) return TGeoShape::Big();
1931 }
1932 // compute distance to get outside this shape
1933 return TGeoArb8::DistFromInside(point, dir, iact, step, safe);
1934}
1935
1936////////////////////////////////////////////////////////////////////////////////
1937/// Compute distance from inside point to surface of the shape.
1938
1939Double_t TGeoGtra::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
1940{
1941 if (iact<3 && safe) {
1942 // compute safe distance
1943 *safe = Safety(point, kTRUE);
1944 if (iact==0) return TGeoShape::Big();
1945 if (iact==1 && step<*safe) return TGeoShape::Big();
1946 }
1947 // compute distance to get outside this shape
1948 return TGeoArb8::DistFromOutside(point, dir, iact, step, safe);
1949}
1950
1951////////////////////////////////////////////////////////////////////////////////
1952/// In case shape has some negative parameters, these has to be computed
1953/// in order to fit the mother
1954
1956{
1957 if (!TestShapeBit(kGeoRunTimeShape)) return 0;
1958 if (mother->IsRunTimeShape()) {
1959 Error("GetMakeRuntimeShape", "invalid mother");
1960 return 0;
1961 }
1962 Double_t dz, h1, bl1, tl1, h2, bl2, tl2;
1963 if (fDz<0) dz=((TGeoTrap*)mother)->GetDz();
1964 else dz=fDz;
1965 if (fH1<0)
1966 h1 = ((TGeoTrap*)mother)->GetH1();
1967 else
1968 h1 = fH1;
1969 if (fH2<0)
1970 h2 = ((TGeoTrap*)mother)->GetH2();
1971 else
1972 h2 = fH2;
1973 if (fBl1<0)
1974 bl1 = ((TGeoTrap*)mother)->GetBl1();
1975 else
1976 bl1 = fBl1;
1977 if (fBl2<0)
1978 bl2 = ((TGeoTrap*)mother)->GetBl2();
1979 else
1980 bl2 = fBl2;
1981 if (fTl1<0)
1982 tl1 = ((TGeoTrap*)mother)->GetTl1();
1983 else
1984 tl1 = fTl1;
1985 if (fTl2<0)
1986 tl2 = ((TGeoTrap*)mother)->GetTl2();
1987 else
1988 tl2 = fTl2;
1989 return (new TGeoGtra(dz, fTheta, fPhi, fTwistAngle ,h1, bl1, tl1, fAlpha1, h2, bl2, tl2, fAlpha2));
1990}
1991
1992////////////////////////////////////////////////////////////////////////////////
1993/// Computes the closest distance from given point to this shape.
1994
1996{
1997 return TGeoArb8::Safety(point,in);
1998}
1999
2000////////////////////////////////////////////////////////////////////////////////
2001/// Save a primitive as a C++ statement(s) on output stream "out".
2002
2003void TGeoGtra::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
2004{
2006 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
2007 out << " dz = " << fDz << ";" << std::endl;
2008 out << " theta = " << fTheta << ";" << std::endl;
2009 out << " phi = " << fPhi << ";" << std::endl;
2010 out << " twist = " << fTwistAngle << ";" << std::endl;
2011 out << " h1 = " << fH1<< ";" << std::endl;
2012 out << " bl1 = " << fBl1<< ";" << std::endl;
2013 out << " tl1 = " << fTl1<< ";" << std::endl;
2014 out << " alpha1 = " << fAlpha1 << ";" << std::endl;
2015 out << " h2 = " << fH2 << ";" << std::endl;
2016 out << " bl2 = " << fBl2<< ";" << std::endl;
2017 out << " tl2 = " << fTl2<< ";" << std::endl;
2018 out << " alpha2 = " << fAlpha2 << ";" << std::endl;
2019 out << " TGeoShape *" << GetPointerName() << " = new TGeoGtra(\"" << GetName() << "\", dz,theta,phi,twist,h1,bl1,tl1,alpha1,h2,bl2,tl2,alpha2);" << std::endl;
2021}
2022
2023////////////////////////////////////////////////////////////////////////////////
2024/// Set all arb8 params in one step.
2025/// - param[0] = dz
2026/// - param[1] = theta
2027/// - param[2] = phi
2028/// - param[3] = h1
2029/// - param[4] = bl1
2030/// - param[5] = tl1
2031/// - param[6] = alpha1
2032/// - param[7] = h2
2033/// - param[8] = bl2
2034/// - param[9] = tl2
2035/// - param[10] = alpha2
2036/// - param[11] = twist
2037
2039{
2041 fTwistAngle = param[11];
2042 Double_t x,y;
2043 Double_t twist = fTwistAngle;
2046 // Coordinates of the center of the bottom face
2047 Double_t xc = -fDz*TMath::Sin(th)*TMath::Cos(ph);
2048 Double_t yc = -fDz*TMath::Sin(th)*TMath::Sin(ph);
2049
2050 Int_t i;
2051
2052 for (i=0; i<4; i++) {
2053 x = fXY[i][0] - xc;
2054 y = fXY[i][1] - yc;
2055 fXY[i][0] = x*TMath::Cos(-0.5*twist*TMath::DegToRad()) + y*TMath::Sin(-0.5*twist*TMath::DegToRad()) + xc;
2056 fXY[i][1] = -x*TMath::Sin(-0.5*twist*TMath::DegToRad()) + y*TMath::Cos(-0.5*twist*TMath::DegToRad()) + yc;
2057 }
2058 // Coordinates of the center of the top face
2059 xc = -xc;
2060 yc = -yc;
2061 for (i=4; i<8; i++) {
2062 x = fXY[i][0] - xc;
2063 y = fXY[i][1] - yc;
2064 fXY[i][0] = x*TMath::Cos(0.5*twist*TMath::DegToRad()) + y*TMath::Sin(0.5*twist*TMath::DegToRad()) + xc;
2065 fXY[i][1] = -x*TMath::Sin(0.5*twist*TMath::DegToRad()) + y*TMath::Cos(0.5*twist*TMath::DegToRad()) + yc;
2066 }
2067 ComputeTwist();
2068 if ((fDz<0) || (fH1<0) || (fBl1<0) || (fTl1<0) ||
2069 (fH2<0) || (fBl2<0) || (fTl2<0)) SetShapeBit(kGeoRunTimeShape);
2070 else TGeoArb8::ComputeBBox();
2071}
2072
2073////////////////////////////////////////////////////////////////////////////////
2074/// Compute distance from array of input points having directions specified by dirs. Store output in dists
2075
2076void TGeoGtra::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2077{
2078 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
2079}
2080
2081////////////////////////////////////////////////////////////////////////////////
2082/// Compute distance from array of input points having directions specified by dirs. Store output in dists
2083
2084void TGeoGtra::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
2085{
2086 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
2087}
2088
2089////////////////////////////////////////////////////////////////////////////////
2090/// Compute safe distance from each of the points in the input array.
2091/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
2092/// Output: Safety values
2093
2094void TGeoGtra::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
2095{
2096 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
2097}
int Int_t
Definition: CPyCppyy.h:43
void Class()
Definition: Class.C:29
#define d(i)
Definition: RSha256.hxx:102
#define b(i)
Definition: RSha256.hxx:100
#define c(i)
Definition: RSha256.hxx:101
static const double x2[5]
static const double x1[5]
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
double Double_t
Definition: RtypesCore.h:59
float Float_t
Definition: RtypesCore.h:57
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:364
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
float xmin
Definition: THbookFile.cxx:95
float ymin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
float ymax
Definition: THbookFile.cxx:95
point * points
Definition: X3DBuffer.c:22
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
Bool_t IsReading() const
Definition: TBuffer.h:86
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
An arbitrary trapezoid with less than 8 vertices standing on two parallel planes perpendicular to Z a...
Definition: TGeoArb8.h:18
virtual void SetVertex(Int_t vnum, Double_t x, Double_t y)
Set values for a given vertex.
Definition: TGeoArb8.cxx:1222
virtual Bool_t GetPointsOnFacet(Int_t, Int_t, Double_t *) const
Fills array with n random points located on the surface of indexed facet.
Definition: TGeoArb8.cxx:932
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this shape along one axis.
Definition: TGeoArb8.cxx:801
Double_t GetClosestEdge(const Double_t *point, Double_t *vert, Int_t &isegment) const
Get index of the edge of the quadrilater represented by vert closest to point.
Definition: TGeoArb8.cxx:371
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get shape range on a given axis.
Definition: TGeoArb8.cxx:811
Double_t fXY[8][2]
[4] tangents of twist angles
Definition: TGeoArb8.h:29
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: TGeoArb8.cxx:1280
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
Computes the closest distance from given point to this shape.
Definition: TGeoArb8.cxx:1023
void CopyTwist(Double_t *twist=nullptr)
Copy twist values from source array.
Definition: TGeoArb8.cxx:214
virtual ~TGeoArb8()
Destructor.
Definition: TGeoArb8.cxx:206
static void GetPlaneNormal(Double_t *p1, Double_t *p2, Double_t *p3, Double_t *norm)
Computes normal to plane defined by P1, P2 and P3.
Definition: TGeoArb8.cxx:904
virtual Bool_t Contains(const Double_t *point) const
Test if point is inside this shape.
Definition: TGeoArb8.cxx:503
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoArb8.cxx:1140
virtual void SetPoints(Double_t *points) const
Creates arb8 mesh points.
Definition: TGeoArb8.cxx:1198
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
Computes distance from outside point to surface of the shape.
Definition: TGeoArb8.cxx:647
void ComputeTwist()
Computes tangents of twist angles (angles between projections on XY plane of corresponding -dz +dz ed...
Definition: TGeoArb8.cxx:271
virtual void Sizeof3D() const
Fill size of this 3-D object.
Definition: TGeoArb8.cxx:1239
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 shape.
Definition: TGeoArb8.cxx:675
Double_t DistToPlane(const Double_t *point, const Double_t *dir, Int_t ipl, Bool_t in) const
Computes distance to plane ipl :
Definition: TGeoArb8.cxx:527
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: TGeoArb8.cxx:1272
static Bool_t InsidePolygon(Double_t x, Double_t y, Double_t *pts)
Finds if a point in XY plane is inside the polygon defines by PTS.
Definition: TGeoArb8.cxx:988
virtual void ComputeBBox()
Computes bounding box for an Arb8 shape.
Definition: TGeoArb8.cxx:245
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: TGeoArb8.cxx:1298
Double_t fDz
Definition: TGeoArb8.h:27
Double_t * fTwist
Definition: TGeoArb8.h:28
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoArb8.cxx:228
Double_t SafetyToFace(const Double_t *point, Int_t iseg, Bool_t in) const
Estimate safety to lateral plane defined by segment iseg in range [0,3] Might be negative: plane seen...
Definition: TGeoArb8.cxx:1101
Bool_t IsTwisted() const
Definition: TGeoArb8.h:75
TGeoArb8()
Default constructor.
Definition: TGeoArb8.cxx:145
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoArb8.cxx:438
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoArb8.cxx:829
Double_t GetTwist(Int_t iseg) const
Get twist for segment I in range [0,3].
Definition: TGeoArb8.cxx:361
static Bool_t IsSamePoint(const Double_t *p1, const Double_t *p2)
Definition: TGeoArb8.h:72
virtual void SetDimensions(Double_t *param)
Set all arb8 params in one step.
Definition: TGeoArb8.cxx:1184
void SetPlaneVertices(Double_t zpl, Double_t *vertices) const
Computes intersection points between plane at zpl and non-horizontal edges.
Definition: TGeoArb8.cxx:1168
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: TGeoArb8.cxx:1262
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: TGeoArb8.cxx:1288
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoArb8.cxx:1008
virtual Int_t GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const
Fills real parameters of a positioned box inside this arb8. Returns 0 if successful.
Definition: TGeoArb8.cxx:848
Box class.
Definition: TGeoBBox.h:18
Double_t fDX
Definition: TGeoBBox.h:21
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
Computes the closest distance from given point to this shape.
Definition: TGeoBBox.cxx:852
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
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoBBox.cxx:790
virtual Double_t GetDX() const
Definition: TGeoBBox.h:74
virtual const Double_t * GetOrigin() const
Definition: TGeoBBox.h:77
virtual Double_t GetDZ() const
Definition: TGeoBBox.h:76
virtual void Sizeof3D() const
Definition: TGeoBBox.cxx:993
virtual Double_t GetDY() const
Definition: TGeoBBox.h:75
Double_t fOrigin[3]
Definition: TGeoBBox.h:24
Double_t fDY
Definition: TGeoBBox.h:22
virtual Bool_t Contains(const Double_t *point) const
Test if point is inside this shape.
Definition: TGeoBBox.cxx:339
Double_t fDZ
Definition: TGeoBBox.h:23
Gtra is a twisted trapezoid.
Definition: TGeoArb8.h:146
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: TGeoArb8.cxx:2094
virtual ~TGeoGtra()
Destructor.
Definition: TGeoArb8.cxx:1917
Double_t fTwistAngle
Definition: TGeoArb8.h:149
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
In case shape has some negative parameters, these has to be computed in order to fit the mother.
Definition: TGeoArb8.cxx:1955
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 inside point to surface of the shape.
Definition: TGeoArb8.cxx:1939
TGeoGtra()
Default ctor.
Definition: TGeoArb8.cxx:1830
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
Computes the closest distance from given point to this shape.
Definition: TGeoArb8.cxx:1995
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: TGeoArb8.cxx:2076
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoArb8.cxx:2003
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 shape.
Definition: TGeoArb8.cxx:1924
virtual void SetDimensions(Double_t *param)
Set all arb8 params in one step.
Definition: TGeoArb8.cxx:2038
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: TGeoArb8.cxx:2084
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Geometrical transformation package.
Definition: TGeoMatrix.h:41
Bool_t IsRotation() const
Definition: TGeoMatrix.h:68
virtual void LocalToMaster(const Double_t *local, Double_t *master) const
convert a point by multiplying its column vector (x, y, z, 1) to matrix inverse
Definition: TGeoMatrix.cxx:339
Node containing an offset.
Definition: TGeoNode.h:184
Base finder class for patterns.
void SetDivIndex(Int_t index)
Base abstract class for all shapes.
Definition: TGeoShape.h:26
static Double_t Big()
Definition: TGeoShape.h:88
static Bool_t IsSegCrossing(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3, Double_t x4, Double_t y4)
Check if segments (A,B) and (C,D) are crossing, where: A(x1,y1), B(x2,y2), C(x3,y3),...
Definition: TGeoShape.cxx:336
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
Definition: TGeoShape.cxx:524
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
Definition: TGeoShape.cxx:326
Bool_t IsRunTimeShape() const
Definition: TGeoShape.h:139
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:699
virtual const char * GetName() const
Get the shape name.
Definition: TGeoShape.cxx:248
@ kGeoClosedShape
Definition: TGeoShape.h:60
@ kGeoSavePrimitive
Definition: TGeoShape.h:65
@ kGeoArb8
Definition: TGeoShape.h:53
@ kGeoRunTimeShape
Definition: TGeoShape.h:41
static Double_t Tolerance()
Definition: TGeoShape.h:91
Bool_t TestShapeBit(UInt_t f) const
Definition: TGeoShape.h:163
TRAP is a general trapezoid, i.e.
Definition: TGeoArb8.h:92
Double_t fTl2
Definition: TGeoArb8.h:103
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: TGeoArb8.cxx:1810
Double_t fBl2
Definition: TGeoArb8.h:102
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this trapezoid shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition: TGeoArb8.cxx:1589
Double_t fAlpha2
Definition: TGeoArb8.h:104
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoArb8.cxx:1733
Double_t fBl1
Definition: TGeoArb8.h:98
Double_t fH1
Definition: TGeoArb8.h:97
virtual TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const
In case shape has some negative parameters, these have to be computed in order to fit the mother.
Definition: TGeoArb8.cxx:1637
Double_t fPhi
Definition: TGeoArb8.h:96
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: TGeoArb8.cxx:1820
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 trapezoid.
Definition: TGeoArb8.cxx:1421
virtual void SetDimensions(Double_t *param)
Set all arb8 params in one step.
Definition: TGeoArb8.cxx:1766
TGeoTrap()
Default ctor.
Definition: TGeoArb8.cxx:1308
Double_t fTl1
Definition: TGeoArb8.h:99
virtual ~TGeoTrap()
Destructor.
Definition: TGeoArb8.cxx:1414
Double_t fH2
Definition: TGeoArb8.h:101
Double_t fTheta
Definition: TGeoArb8.h:95
virtual Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const
Computes the closest distance from given point to this shape.
Definition: TGeoArb8.cxx:1672
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: TGeoArb8.cxx:1802
Double_t fAlpha1
Definition: TGeoArb8.h:100
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 trapezoid.
Definition: TGeoArb8.cxx:1474
Volume families.
Definition: TGeoVolume.h:254
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:49
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
TGeoMedium * GetMedium() const
Definition: TGeoVolume.h:173
void SetFinder(TGeoPatternFinder *finder)
Definition: TGeoVolume.h:231
Int_t GetNdaughters() const
Definition: TGeoVolume.h:349
TObjArray * GetNodes()
Definition: TGeoVolume.h:167
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition: TObject.cxx:921
Basic string class.
Definition: TString.h:136
const char * Data() const
Definition: TString.h:369
TPaveText * pt
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
TH1F * h1
Definition: legend1.C:5
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
static const std::string name("name")
static constexpr double s
Long64_t LocMin(Long64_t n, const T *a)
Return index of array with the minimum element.
Definition: TMath.h:972
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:212
T1 Sign(T1 a, T2 b)
Definition: TMathBase.h:165
Long64_t LocMax(Long64_t n, const T *a)
Return index of array with the maximum element.
Definition: TMath.h:1000
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:96
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:81
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
Double_t Cos(Double_t)
Definition: TMath.h:643
Double_t Sin(Double_t)
Definition: TMath.h:639
Double_t Tan(Double_t)
Definition: TMath.h:647
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
auto * a
Definition: textangle.C:12
#define snext(osub1, osub2)
Definition: triangle.c:1168