Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoTrd1.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Andrei Gheata 24/10/01
3// TGeoTrd1::Contains() and DistFromInside() implemented by Mihaela Gheata
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13/** \class TGeoTrd1
14\ingroup Trapezoids
15
16A trapezoid with only X varying with Z. It is defined by the
17half-length in Z, the half-length in X at the lowest and highest Z
18planes and the half-length in Y:
19
20~~~ {.cpp}
21TGeoTrd1(Double_t dx1,Double_t dx2,Double_t dy,Double_t dz);
22~~~
23
24Begin_Macro
25{
26 TCanvas *c = new TCanvas("c", "c",0,0,600,600);
27 new TGeoManager("trd1", "poza8");
28 TGeoMaterial *mat = new TGeoMaterial("Al", 26.98,13,2.7);
29 TGeoMedium *med = new TGeoMedium("MED",1,mat);
30 TGeoVolume *top = gGeoManager->MakeBox("TOP",med,100,100,100);
31 gGeoManager->SetTopVolume(top);
32 TGeoVolume *vol = gGeoManager->MakeTrd1("Trd1",med, 10,20,30,40);
33 vol->SetLineWidth(2);
34 top->AddNode(vol,1);
35 gGeoManager->CloseGeometry();
36 gGeoManager->SetNsegments(80);
37 top->Draw();
38 TView *view = gPad->GetView();
39 if (view) view->ShowAxis();
40}
41End_Macro
42*/
43
44#include <iostream>
45
46#include "TGeoManager.h"
47#include "TGeoMatrix.h"
48#include "TGeoVolume.h"
49#include "TGeoTrd1.h"
50#include "TMath.h"
51
52
53////////////////////////////////////////////////////////////////////////////////
54/// dummy ctor
55
57{
58 fDz = fDx1 = fDx2 = fDy = 0;
60}
61
62////////////////////////////////////////////////////////////////////////////////
63/// constructor.
64
66{
68 fDx1 = dx1;
69 fDx2 = dx2;
70 fDy = dy;
71 fDz = dz;
72 if ((dx1 < 0) || (dx2 < 0) || (dy < 0) || (dz < 0)) {
74 printf("trd1 : dx1=%f, dx2=%f, dy=%f, dz=%f\n", dx1, dx2, dy, dz);
75 } else
77}
78
79////////////////////////////////////////////////////////////////////////////////
80/// constructor.
81
83{
85 fDx1 = dx1;
86 fDx2 = dx2;
87 fDy = dy;
88 fDz = dz;
89 if ((dx1 < 0) || (dx2 < 0) || (dy < 0) || (dz < 0)) {
91 printf("trd1 : dx1=%f, dx2=%f, dy=%f, dz=%f\n", dx1, dx2, dy, dz);
92 } else
94}
95
96////////////////////////////////////////////////////////////////////////////////
97/// ctor with an array of parameters
98/// - param[0] = dx1
99/// - param[1] = dx2
100/// - param[2] = dy
101/// - param[3] = dz
102
104{
106 SetDimensions(param);
107 if ((fDx1 < 0) || (fDx2 < 0) || (fDy <= 0) || (fDz <= 0))
109 else
110 ComputeBBox();
111}
112
113////////////////////////////////////////////////////////////////////////////////
114/// destructor
115
117
118////////////////////////////////////////////////////////////////////////////////
119/// Computes capacity of the shape in [length^3]
120
122{
123 Double_t capacity = 4. * (fDx1 + fDx2) * fDy * fDz;
124 return capacity;
125}
126
127////////////////////////////////////////////////////////////////////////////////
128/// compute bounding box for a trd1
129
131{
133 fDY = fDy;
134 fDZ = fDz;
135 memset(fOrigin, 0, 3 * sizeof(Double_t));
136}
137
138////////////////////////////////////////////////////////////////////////////////
139/// Compute normal to closest surface from POINT.
140
141void TGeoTrd1::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const
142{
144 //--- Compute safety first
145 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
146 Double_t calf = 1. / TMath::Sqrt(1.0 + fx * fx);
147 // check Z facettes
148 safe = safemin = TMath::Abs(fDz - TMath::Abs(point[2]));
149 norm[0] = norm[1] = 0;
150 norm[2] = (dir[2] >= 0) ? 1 : -1;
151 if (safe < 1E-6)
152 return;
153 // check X facettes
154 Double_t distx = 0.5 * (fDx1 + fDx2) - fx * point[2];
155 if (distx >= 0) {
156 safe = TMath::Abs(distx - TMath::Abs(point[0])) * calf;
157 if (safe < safemin) {
158 safemin = safe;
159 norm[0] = (point[0] > 0) ? calf : (-calf);
160 norm[1] = 0;
161 norm[2] = calf * fx;
162 Double_t dot = norm[0] * dir[0] + norm[1] * dir[1] + norm[2] * dir[2];
163 if (dot < 0) {
164 norm[0] = -norm[0];
165 norm[2] = -norm[2];
166 }
167 if (safe < 1E-6)
168 return;
169 }
170 }
171 // check Y facettes
172 safe = TMath::Abs(fDy - TMath::Abs(point[1]));
173 if (safe < safemin) {
174 norm[0] = norm[2] = 0;
175 norm[1] = (dir[1] >= 0) ? 1 : -1;
176 }
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// test if point is inside this shape
181/// check Z range
182
184{
185 if (TMath::Abs(point[2]) > fDz)
186 return kFALSE;
187 // then y
188 if (TMath::Abs(point[1]) > fDy)
189 return kFALSE;
190 // then x
191 Double_t dx = 0.5 * (fDx2 * (point[2] + fDz) + fDx1 * (fDz - point[2])) / fDz;
192 if (TMath::Abs(point[0]) > dx)
193 return kFALSE;
194 return kTRUE;
195}
196
197////////////////////////////////////////////////////////////////////////////////
198/// Compute distance from inside point to surface of the trd1
199/// Boundary safe algorithm.
200
203{
204 if (iact < 3 && safe) {
205 // compute safe distance
206 *safe = Safety(point, kTRUE);
207 if (iact == 0)
208 return TGeoShape::Big();
209 if (iact == 1 && step < *safe)
210 return TGeoShape::Big();
211 }
212
213 //--- Compute safety first
214 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
215 Double_t cn;
216 Double_t distx = 0.5 * (fDx1 + fDx2) - fx * point[2];
217 //--- Compute distance to this shape
218 // first check if Z facettes are crossed
219 Double_t dist[3];
220 for (Int_t i = 0; i < 3; i++)
221 dist[i] = TGeoShape::Big();
222 if (dir[2] < 0) {
223 dist[0] = -(point[2] + fDz) / dir[2];
224 } else if (dir[2] > 0) {
225 dist[0] = (fDz - point[2]) / dir[2];
226 }
227 if (dist[0] <= 0)
228 return 0.0;
229 // now check X facettes
230 cn = -dir[0] + fx * dir[2];
231 if (cn > 0) {
232 dist[1] = point[0] + distx;
233 if (dist[1] <= 0)
234 return 0.0;
235 dist[1] /= cn;
236 }
237 cn = dir[0] + fx * dir[2];
238 if (cn > 0) {
239 Double_t s = distx - point[0];
240 if (s <= 0)
241 return 0.0;
242 s /= cn;
243 if (s < dist[1])
244 dist[1] = s;
245 }
246 // now check Y facettes
247 if (dir[1] < 0) {
248 dist[2] = -(point[1] + fDy) / dir[1];
249 } else if (dir[1] > 0) {
250 dist[2] = (fDy - point[1]) / dir[1];
251 }
252 if (dist[2] <= 0)
253 return 0.0;
254 return dist[TMath::LocMin(3, dist)];
255}
256
257////////////////////////////////////////////////////////////////////////////////
258/// get the most visible corner from outside point and the normals
259
261{
262 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
263 Double_t calf = 1. / TMath::Sqrt(1.0 + fx * fx);
264 Double_t salf = calf * fx;
265 // check visibility of X faces
266 Double_t distx = 0.5 * (fDx1 + fDx2) - fx * point[2];
267 memset(normals, 0, 9 * sizeof(Double_t));
268 TGeoTrd1 *trd1 = (TGeoTrd1 *)this;
269 if (point[0] > distx) {
270 // hi x face visible
271 trd1->SetShapeBit(kGeoVisX);
272 normals[0] = calf;
273 normals[2] = salf;
274 } else {
275 trd1->SetShapeBit(kGeoVisX, kFALSE);
276 normals[0] = -calf;
277 normals[2] = salf;
278 }
279 if (point[1] > fDy) {
280 // hi y face visible
281 trd1->SetShapeBit(kGeoVisY);
282 normals[4] = 1;
283 } else {
284 trd1->SetShapeBit(kGeoVisY, kFALSE);
285 normals[4] = -1;
286 }
287 if (point[2] > fDz) {
288 // hi z face visible
289 trd1->SetShapeBit(kGeoVisZ);
290 normals[8] = 1;
291 } else {
292 trd1->SetShapeBit(kGeoVisZ, kFALSE);
293 normals[8] = -1;
294 }
296}
297
298////////////////////////////////////////////////////////////////////////////////
299/// get the opposite corner of the intersected face
300
302{
303 TGeoTrd1 *trd1 = (TGeoTrd1 *)this;
304 if (inorm != 0) {
305 // change x face
306 trd1->SetShapeBit(kGeoVisX, !TestShapeBit(kGeoVisX));
307 normals[0] = -normals[0];
308 }
309 if (inorm != 1) {
310 // change y face
311 trd1->SetShapeBit(kGeoVisY, !TestShapeBit(kGeoVisY));
312 normals[4] = -normals[4];
313 }
314 if (inorm != 2) {
315 // hi z face visible
316 trd1->SetShapeBit(kGeoVisZ, !TestShapeBit(kGeoVisZ));
317 normals[8] = -normals[8];
318 }
320}
321
322////////////////////////////////////////////////////////////////////////////////
323/// Compute distance from outside point to surface of the trd1
324/// Boundary safe algorithm
325
328{
329 if (iact < 3 && safe) {
330 // compute safe distance
331 *safe = Safety(point, kFALSE);
332 if (iact == 0)
333 return TGeoShape::Big();
334 if (iact == 1 && step < *safe)
335 return TGeoShape::Big();
336 }
337 // find a visible face
339 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
340 Double_t cn;
341 Double_t distx = 0.5 * (fDx1 + fDx2) - fx * point[2];
342 Bool_t in = kTRUE;
343 Double_t safx = distx - TMath::Abs(point[0]);
344 Double_t safy = fDy - TMath::Abs(point[1]);
345 Double_t safz = fDz - TMath::Abs(point[2]);
346
347 //--- Compute distance to this shape
348 // first check if Z facettes are crossed
349 if (point[2] <= -fDz) {
350 if (dir[2] <= 0)
351 return TGeoShape::Big();
352 in = kFALSE;
353 Double_t snxt = -(fDz + point[2]) / dir[2];
354 // find extrapolated X and Y
355 xnew = point[0] + snxt * dir[0];
356 if (TMath::Abs(xnew) <= fDx1) {
357 ynew = point[1] + snxt * dir[1];
358 if (TMath::Abs(ynew) <= fDy)
359 return snxt;
360 }
361 } else if (point[2] >= fDz) {
362 if (dir[2] >= 0)
363 return TGeoShape::Big();
364 in = kFALSE;
365 Double_t snxt = (fDz - point[2]) / dir[2];
366 // find extrapolated X and Y
367 xnew = point[0] + snxt * dir[0];
368 if (TMath::Abs(xnew) <= fDx2) {
369 ynew = point[1] + snxt * dir[1];
370 if (TMath::Abs(ynew) <= fDy)
371 return snxt;
372 }
373 }
374 // check if X facettes are crossed
375 if (point[0] <= -distx) {
376 cn = -dir[0] + fx * dir[2];
377 if (cn >= 0)
378 return TGeoShape::Big();
379 in = kFALSE;
380 Double_t snxt = (point[0] + distx) / cn;
381 // find extrapolated Y and Z
382 ynew = point[1] + snxt * dir[1];
383 if (TMath::Abs(ynew) <= fDy) {
384 znew = point[2] + snxt * dir[2];
385 if (TMath::Abs(znew) <= fDz)
386 return snxt;
387 }
388 }
389 if (point[0] >= distx) {
390 cn = dir[0] + fx * dir[2];
391 if (cn >= 0)
392 return TGeoShape::Big();
393 in = kFALSE;
394 Double_t snxt = (distx - point[0]) / cn;
395 // find extrapolated Y and Z
396 ynew = point[1] + snxt * dir[1];
397 if (TMath::Abs(ynew) < fDy) {
398 znew = point[2] + snxt * dir[2];
399 if (TMath::Abs(znew) < fDz)
400 return snxt;
401 }
402 }
403 // finally check Y facettes
404 if (point[1] <= -fDy) {
405 cn = -dir[1];
406 if (cn >= 0)
407 return TGeoShape::Big();
408 in = kFALSE;
409 Double_t snxt = (point[1] + fDy) / cn;
410 // find extrapolated X and Z
411 znew = point[2] + snxt * dir[2];
412 if (TMath::Abs(znew) < fDz) {
413 xnew = point[0] + snxt * dir[0];
414 Double_t dx = 0.5 * (fDx1 + fDx2) - fx * znew;
415 if (TMath::Abs(xnew) < dx)
416 return snxt;
417 }
418 } else if (point[1] >= fDy) {
419 cn = dir[1];
420 if (cn >= 0)
421 return TGeoShape::Big();
422 in = kFALSE;
423 Double_t snxt = (fDy - point[1]) / cn;
424 // find extrapolated X and Z
425 znew = point[2] + snxt * dir[2];
426 if (TMath::Abs(znew) < fDz) {
427 xnew = point[0] + snxt * dir[0];
428 Double_t dx = 0.5 * (fDx1 + fDx2) - fx * znew;
429 if (TMath::Abs(xnew) < dx)
430 return snxt;
431 }
432 }
433 if (!in)
434 return TGeoShape::Big();
435 // Point actually inside
436 if (safz < safx && safz < safy) {
437 if (point[2] * dir[2] >= 0)
438 return TGeoShape::Big();
439 return 0.0;
440 }
441 if (safy < safx) {
442 if (point[1] * dir[1] >= 0)
443 return TGeoShape::Big();
444 return 0.0;
445 }
446 cn = TMath::Sign(1.0, point[0]) * dir[0] + fx * dir[2];
447 if (cn >= 0)
448 return TGeoShape::Big();
449 return 0.0;
450}
451
452////////////////////////////////////////////////////////////////////////////////
453/// Divide this trd1 shape belonging to volume "voldiv" into ndiv volumes
454/// called divname, from start position with the given step. Returns pointer
455/// to created division cell volume in case of Y divisions. For Z divisions just
456/// return the pointer to the volume to be divided. In case a wrong
457/// division axis is supplied, returns pointer to volume that was divided.
458
461{
462 TGeoShape *shape; //--- shape to be created
463 TGeoVolume *vol; //--- division volume to be created
464 TGeoVolumeMulti *vmulti; //--- generic divided volume
465 TGeoPatternFinder *finder; //--- finder to be attached
466 TString opt = ""; //--- option to be attached
467 Double_t zmin, zmax, dx1n, dx2n;
468 Int_t id;
469 Double_t end = start + ndiv * step;
470 switch (iaxis) {
471 case 1: Warning("Divide", "dividing a Trd1 on X not implemented"); return nullptr;
472 case 2:
473 finder = new TGeoPatternY(voldiv, ndiv, start, end);
474 voldiv->SetFinder(finder);
475 finder->SetDivIndex(voldiv->GetNdaughters());
476 shape = new TGeoTrd1(fDx1, fDx2, step / 2, fDz);
477 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
479 vmulti->AddVolume(vol);
480 opt = "Y";
481 for (id = 0; id < ndiv; id++) {
482 voldiv->AddNodeOffset(vol, id, start + step / 2 + id * step, opt.Data());
483 ((TGeoNodeOffset *)voldiv->GetNodes()->At(voldiv->GetNdaughters() - 1))->SetFinder(finder);
484 }
485 return vmulti;
486 case 3:
487 finder = new TGeoPatternZ(voldiv, ndiv, start, end);
488 voldiv->SetFinder(finder);
489 finder->SetDivIndex(voldiv->GetNdaughters());
491 for (id = 0; id < ndiv; id++) {
492 zmin = start + id * step;
493 zmax = start + (id + 1) * step;
494 dx1n = 0.5 * (fDx1 * (fDz - zmin) + fDx2 * (fDz + zmin)) / fDz;
495 dx2n = 0.5 * (fDx1 * (fDz - zmax) + fDx2 * (fDz + zmax)) / fDz;
496 shape = new TGeoTrd1(dx1n, dx2n, fDy, step / 2.);
497 vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
498 vmulti->AddVolume(vol);
499 opt = "Z";
500 voldiv->AddNodeOffset(vol, id, start + step / 2 + id * step, opt.Data());
501 ((TGeoNodeOffset *)voldiv->GetNodes()->At(voldiv->GetNdaughters() - 1))->SetFinder(finder);
502 }
503 return vmulti;
504 default: Error("Divide", "Wrong axis type for division"); return nullptr;
505 }
506}
507
508////////////////////////////////////////////////////////////////////////////////
509/// Get range of shape for a given axis.
510
512{
513 xlo = 0;
514 xhi = 0;
515 Double_t dx = 0;
516 switch (iaxis) {
517 case 2:
518 xlo = -fDy;
519 xhi = fDy;
520 dx = xhi - xlo;
521 return dx;
522 case 3:
523 xlo = -fDz;
524 xhi = fDz;
525 dx = xhi - xlo;
526 return dx;
527 }
528 return dx;
529}
530
531////////////////////////////////////////////////////////////////////////////////
532/// Fill vector param[4] with the bounding cylinder parameters. The order
533/// is the following : Rmin, Rmax, Phi1, Phi2
534
539
540////////////////////////////////////////////////////////////////////////////////
541/// Fills real parameters of a positioned box inside this. Returns 0 if successful.
542
544{
545 dx = dy = dz = 0;
546 if (mat->IsRotation()) {
547 Error("GetFittingBox", "cannot handle parametrized rotated volumes");
548 return 1; // ### rotation not accepted ###
549 }
550 //--> translate the origin of the parametrized box to the frame of this box.
551 Double_t origin[3];
552 mat->LocalToMaster(parambox->GetOrigin(), origin);
553 if (!Contains(origin)) {
554 Error("GetFittingBox", "wrong matrix - parametrized box is outside this");
555 return 1; // ### wrong matrix ###
556 }
557 //--> now we have to get the valid range for all parametrized axis
558 Double_t dd[3];
559 dd[0] = parambox->GetDX();
560 dd[1] = parambox->GetDY();
561 dd[2] = parambox->GetDZ();
562 //-> check if Z range is fixed
563 if (dd[2] < 0) {
564 dd[2] = TMath::Min(origin[2] + fDz, fDz - origin[2]);
565 if (dd[2] < 0) {
566 Error("GetFittingBox", "wrong matrix");
567 return 1;
568 }
569 }
570 //-> check if Y range is fixed
571 if (dd[1] < 0) {
572 dd[1] = TMath::Min(origin[1] + fDy, fDy - origin[1]);
573 if (dd[1] < 0) {
574 Error("GetFittingBox", "wrong matrix");
575 return 1;
576 }
577 }
578 if (dd[0] >= 0) {
579 dx = dd[0];
580 dy = dd[1];
581 dz = dd[2];
582 return 0;
583 }
584 //-> check now range at Z = origin[2] +/- dd[2]
585 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
586 Double_t dx0 = 0.5 * (fDx1 + fDx2);
587 Double_t z = origin[2] - dd[2];
588 dd[0] = dx0 - fx * z - origin[0];
589 z = origin[2] + dd[2];
590 dd[0] = TMath::Min(dd[0], dx0 - fx * z - origin[0]);
591 if (dd[0] < 0) {
592 Error("GetFittingBox", "wrong matrix");
593 return 1;
594 }
595 dx = dd[0];
596 dy = dd[1];
597 dz = dd[2];
598 return 0;
599}
600
601////////////////////////////////////////////////////////////////////////////////
602/// in case shape has some negative parameters, these has to be computed
603/// in order to fit the mother
604
606{
608 return nullptr;
609 if (!mother->TestShapeBit(kGeoTrd1)) {
610 Error("GetMakeRuntimeShape", "invalid mother");
611 return nullptr;
612 }
613 Double_t dx1, dx2, dy, dz;
614 if (fDx1 < 0)
615 dx1 = ((TGeoTrd1 *)mother)->GetDx1();
616 else
617 dx1 = fDx1;
618 if (fDx2 < 0)
619 dx2 = ((TGeoTrd1 *)mother)->GetDx2();
620 else
621 dx2 = fDx2;
622 if (fDy < 0)
623 dy = ((TGeoTrd1 *)mother)->GetDy();
624 else
625 dy = fDy;
626 if (fDz < 0)
627 dz = ((TGeoTrd1 *)mother)->GetDz();
628 else
629 dz = fDz;
630
631 return (new TGeoTrd1(dx1, dx2, dy, dz));
632}
633
634////////////////////////////////////////////////////////////////////////////////
635/// print shape parameters
636
638{
639 printf("*** Shape %s: TGeoTrd1 ***\n", GetName());
640 printf(" dx1 = %11.5f\n", fDx1);
641 printf(" dx2 = %11.5f\n", fDx2);
642 printf(" dy = %11.5f\n", fDy);
643 printf(" dz = %11.5f\n", fDz);
644 printf(" Bounding box:\n");
646}
647
648////////////////////////////////////////////////////////////////////////////////
649/// computes the closest distance from given point to this shape, according
650/// to option. The matching point on the shape is stored in spoint.
651
653{
654 Double_t saf[3];
655 //--- Compute safety first
656 // check Z facettes
657 saf[0] = fDz - TMath::Abs(point[2]);
658 Double_t fx = 0.5 * (fDx1 - fDx2) / fDz;
659 Double_t calf = 1. / TMath::Sqrt(1.0 + fx * fx);
660 // check X facettes
661 Double_t distx = 0.5 * (fDx1 + fDx2) - fx * point[2];
662 if (distx < 0)
663 saf[1] = TGeoShape::Big();
664 else
665 saf[1] = (distx - TMath::Abs(point[0])) * calf;
666 // check Y facettes
667 saf[2] = fDy - TMath::Abs(point[1]);
668 if (in)
669 return saf[TMath::LocMin(3, saf)];
670 for (Int_t i = 0; i < 3; i++)
671 saf[i] = -saf[i];
672 return saf[TMath::LocMax(3, saf)];
673}
674
675////////////////////////////////////////////////////////////////////////////////
676/// Save a primitive as a C++ statement(s) on output stream "out".
677
678void TGeoTrd1::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
679{
681 return;
682 out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
683 out << " dx1 = " << fDx1 << ";" << std::endl;
684 out << " dx2 = " << fDx2 << ";" << std::endl;
685 out << " dy = " << fDy << ";" << std::endl;
686 out << " dz = " << fDZ << ";" << std::endl;
687 out << " TGeoShape *" << GetPointerName() << " = new TGeoTrd1(\"" << GetName() << "\", dx1,dx2,dy,dz);"
688 << std::endl;
690}
691
692////////////////////////////////////////////////////////////////////////////////
693/// set trd1 params in one step :
694
696{
697 fDx1 = param[0];
698 fDx2 = param[1];
699 fDy = param[2];
700 fDz = param[3];
701 ComputeBBox();
702}
703
704////////////////////////////////////////////////////////////////////////////////
705/// set vertex of a corner according to visibility flags
706
708{
709 if (TestShapeBit(kGeoVisX)) {
710 if (TestShapeBit(kGeoVisZ)) {
711 vertex[0] = fDx2;
712 vertex[2] = fDz;
713 vertex[1] = (TestShapeBit(kGeoVisY)) ? fDy : -fDy;
714 } else {
715 vertex[0] = fDx1;
716 vertex[2] = -fDz;
717 vertex[1] = (TestShapeBit(kGeoVisY)) ? fDy : -fDy;
718 }
719 } else {
720 if (TestShapeBit(kGeoVisZ)) {
721 vertex[0] = -fDx2;
722 vertex[2] = fDz;
723 vertex[1] = (TestShapeBit(kGeoVisY)) ? fDy : -fDy;
724 } else {
725 vertex[0] = -fDx1;
726 vertex[2] = -fDz;
727 vertex[1] = (TestShapeBit(kGeoVisY)) ? fDy : -fDy;
728 }
729 }
730}
731
732////////////////////////////////////////////////////////////////////////////////
733/// create arb8 mesh points
734
736{
737 if (!points)
738 return;
739 points[0] = -fDx1;
740 points[1] = -fDy;
741 points[2] = -fDz;
742 points[3] = -fDx1;
743 points[4] = fDy;
744 points[5] = -fDz;
745 points[6] = fDx1;
746 points[7] = fDy;
747 points[8] = -fDz;
748 points[9] = fDx1;
749 points[10] = -fDy;
750 points[11] = -fDz;
751 points[12] = -fDx2;
752 points[13] = -fDy;
753 points[14] = fDz;
754 points[15] = -fDx2;
755 points[16] = fDy;
756 points[17] = fDz;
757 points[18] = fDx2;
758 points[19] = fDy;
759 points[20] = fDz;
760 points[21] = fDx2;
761 points[22] = -fDy;
762 points[23] = fDz;
763}
764
765////////////////////////////////////////////////////////////////////////////////
766/// create arb8 mesh points
767
769{
770 if (!points)
771 return;
772 points[0] = -fDx1;
773 points[1] = -fDy;
774 points[2] = -fDz;
775 points[3] = -fDx1;
776 points[4] = fDy;
777 points[5] = -fDz;
778 points[6] = fDx1;
779 points[7] = fDy;
780 points[8] = -fDz;
781 points[9] = fDx1;
782 points[10] = -fDy;
783 points[11] = -fDz;
784 points[12] = -fDx2;
785 points[13] = -fDy;
786 points[14] = fDz;
787 points[15] = -fDx2;
788 points[16] = fDy;
789 points[17] = fDz;
790 points[18] = fDx2;
791 points[19] = fDy;
792 points[20] = fDz;
793 points[21] = fDx2;
794 points[22] = -fDy;
795 points[23] = fDz;
796}
797
798////////////////////////////////////////////////////////////////////////////////
799/// fill size of this 3-D object
800
802{
804}
805
806////////////////////////////////////////////////////////////////////////////////
807/// Check the inside status for each of the points in the array.
808/// Input: Array of point coordinates + vector size
809/// Output: Array of Booleans for the inside of each point
810
812{
813 for (Int_t i = 0; i < vecsize; i++)
814 inside[i] = Contains(&points[3 * i]);
815}
816
817////////////////////////////////////////////////////////////////////////////////
818/// Compute the normal for an array o points so that norm.dot.dir is positive
819/// Input: Arrays of point coordinates and directions + vector size
820/// Output: Array of normal directions
821
823{
824 for (Int_t i = 0; i < vecsize; i++)
825 ComputeNormal(&points[3 * i], &dirs[3 * i], &norms[3 * i]);
826}
827
828////////////////////////////////////////////////////////////////////////////////
829/// Compute distance from array of input points having directions specified by dirs. Store output in dists
830
832 Double_t *step) const
833{
834 for (Int_t i = 0; i < vecsize; i++)
835 dists[i] = DistFromInside(&points[3 * i], &dirs[3 * i], 3, step[i]);
836}
837
838////////////////////////////////////////////////////////////////////////////////
839/// Compute distance from array of input points having directions specified by dirs. Store output in dists
840
842 Double_t *step) const
843{
844 for (Int_t i = 0; i < vecsize; i++)
845 dists[i] = DistFromOutside(&points[3 * i], &dirs[3 * i], 3, step[i]);
846}
847
848////////////////////////////////////////////////////////////////////////////////
849/// Compute safe distance from each of the points in the input array.
850/// Input: Array of point coordinates, array of statuses for these points, size of the arrays
851/// Output: Safety values
852
853void TGeoTrd1::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
854{
855 for (Int_t i = 0; i < vecsize; i++)
856 safe[i] = Safety(&points[3 * i], inside[i]);
857}
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
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
Box class.
Definition TGeoBBox.h:17
Double_t fDX
Definition TGeoBBox.h:20
void GetBoundingCylinder(Double_t *param) const override
Fill vector param[4] with the bounding cylinder parameters.
Definition TGeoBBox.cxx:600
Double_t fOrigin[3]
Definition TGeoBBox.h:23
void InspectShape() const override
Prints shape parameters.
Definition TGeoBBox.cxx:810
void Sizeof3D() const override
Double_t fDY
Definition TGeoBBox.h:21
Double_t fDZ
Definition TGeoBBox.h:22
TGeoVolumeMulti * MakeVolumeMulti(const char *name, TGeoMedium *medium)
Make a TGeoVolumeMulti handling a list of volumes.
Geometrical transformation package.
Definition TGeoMatrix.h:38
Node containing an offset.
Definition TGeoNode.h:184
base finder class for patterns. A pattern is specifying a division type
a Y axis divison pattern
a Z axis divison pattern
Base abstract class for all shapes.
Definition TGeoShape.h:25
static Double_t Big()
Definition TGeoShape.h:94
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
const char * GetPointerName() const
Provide a pointer name containing uid.
const char * GetName() const override
Get the shape name.
@ kGeoSavePrimitive
Definition TGeoShape.h:64
@ kGeoRunTimeShape
Definition TGeoShape.h:40
Bool_t TestShapeBit(UInt_t f) const
Definition TGeoShape.h:175
A trapezoid with only X varying with Z.
Definition TGeoTrd1.h:17
void ComputeBBox() override
compute bounding box for a trd1
Definition TGeoTrd1.cxx:130
void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize) override
Compute the normal for an array o points so that norm.dot.dir is positive Input: Arrays of point coor...
Definition TGeoTrd1.cxx:822
void GetBoundingCylinder(Double_t *param) const override
Fill vector param[4] with the bounding cylinder parameters.
Definition TGeoTrd1.cxx:535
void InspectShape() const override
print shape parameters
Definition TGeoTrd1.cxx:637
TGeoShape * GetMakeRuntimeShape(TGeoShape *mother, TGeoMatrix *mat) const override
in case shape has some negative parameters, these has to be computed in order to fit the mother
Definition TGeoTrd1.cxx:605
Bool_t Contains(const Double_t *point) const override
test if point is inside this shape check Z range
Definition TGeoTrd1.cxx:183
void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
Compute distance from array of input points having directions specified by dirs. Store output in dist...
Definition TGeoTrd1.cxx:831
void SetDimensions(Double_t *param) override
set trd1 params in one step :
Definition TGeoTrd1.cxx:695
void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const override
Check the inside status for each of the points in the array.
Definition TGeoTrd1.cxx:811
void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
Compute distance from array of input points having directions specified by dirs. Store output in dist...
Definition TGeoTrd1.cxx:841
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const override
computes the closest distance from given point to this shape, according to option.
Definition TGeoTrd1.cxx:652
void SetVertex(Double_t *vertex) const
set vertex of a corner according to visibility flags
Definition TGeoTrd1.cxx:707
Int_t GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const override
Fills real parameters of a positioned box inside this. Returns 0 if successful.
Definition TGeoTrd1.cxx:543
void GetOppositeCorner(const Double_t *point, Int_t inorm, Double_t *vertex, Double_t *normals) const
get the opposite corner of the intersected face
Definition TGeoTrd1.cxx:301
void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const override
Compute safe distance from each of the points in the input array.
Definition TGeoTrd1.cxx:853
void GetVisibleCorner(const Double_t *point, Double_t *vertex, Double_t *normals) const
get the most visible corner from outside point and the normals
Definition TGeoTrd1.cxx:260
~TGeoTrd1() override
destructor
Definition TGeoTrd1.cxx:116
Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const override
Get range of shape for a given axis.
Definition TGeoTrd1.cxx:511
void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) const override
Compute normal to closest surface from POINT.
Definition TGeoTrd1.cxx:141
Double_t fDz
Definition TGeoTrd1.h:23
void Sizeof3D() const override
fill size of this 3-D object
Definition TGeoTrd1.cxx:801
TGeoTrd1()
dummy ctor
Definition TGeoTrd1.cxx:56
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
Definition TGeoTrd1.cxx:678
Double_t fDx1
Definition TGeoTrd1.h:20
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
Compute distance from inside point to surface of the trd1 Boundary safe algorithm.
Definition TGeoTrd1.cxx:202
void SetPoints(Double_t *points) const override
create arb8 mesh points
Definition TGeoTrd1.cxx:735
TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step) override
Divide this trd1 shape belonging to volume "voldiv" into ndiv volumes called divname,...
Definition TGeoTrd1.cxx:460
Double_t Capacity() const override
Computes capacity of the shape in [length^3].
Definition TGeoTrd1.cxx:121
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
Compute distance from outside point to surface of the trd1 Boundary safe algorithm.
Definition TGeoTrd1.cxx:327
Double_t fDy
Definition TGeoTrd1.h:22
Double_t fDx2
Definition TGeoTrd1.h:21
Volume families.
Definition TGeoVolume.h:266
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition TGeoVolume.h:43
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:226
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
Long64_t LocMin(Long64_t n, const T *a)
Returns index of array with the minimum element.
Definition TMath.h:993
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
T1 Sign(T1 a, T2 b)
Returns a value with the magnitude of a and the sign of b.
Definition TMathBase.h:176
Long64_t LocMax(Long64_t n, const T *a)
Returns index of array with the maximum element.
Definition TMath.h:1099
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:673
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:199
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124