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