Logo ROOT  
Reference Guide
TGeoPara.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 31/01/02
3 // TGeoPara::Contains() implemented by Mihaela Gheata
4 
5 /*************************************************************************
6  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7  * All rights reserved. *
8  * *
9  * For the licensing terms see $ROOTSYS/LICENSE. *
10  * For the list of contributors see $ROOTSYS/README/CREDITS. *
11  *************************************************************************/
12 
13 /** \class TGeoPara
14 \ingroup Geometry_classes
15 
16 Parallelepiped class. It has 6 parameters :
17 
18  - dx, dy, dz - half lengths in X, Y, Z
19  - alpha - angle w.r.t the Y axis from center of low Y edge to
20  center of high Y edge [deg]
21  - theta, phi - polar and azimuthal angles of the segment between
22  low and high Z surfaces [deg]
23 
24 Begin_Macro(source)
25 {
26  TCanvas *c = new TCanvas("c", "c",0,0,600,600);
27  new TGeoManager("para", "poza1");
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->MakePara("PARA",med, 20,30,40,30,15,30);
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  view->ShowAxis();
40 }
41 End_Macro
42 */
43 
44 #include <iostream>
45 
46 #include "TGeoManager.h"
47 #include "TGeoMatrix.h"
48 #include "TGeoVolume.h"
49 #include "TGeoPara.h"
50 #include "TMath.h"
51 
53 
54 ////////////////////////////////////////////////////////////////////////////////
55 /// Default constructor
56 
58 {
60  fX = fY = fZ = 0;
61  fAlpha = 0;
62  fTheta = 0;
63  fPhi = 0;
64  fTxy = 0;
65  fTxz = 0;
66  fTyz = 0;
67 }
68 
69 ////////////////////////////////////////////////////////////////////////////////
70 /// Default constructor specifying minimum and maximum radius
71 
73  Double_t theta, Double_t phi)
74  :TGeoBBox(0, 0, 0)
75 {
77  fX = dx;
78  fY = dy;
79  fZ = dz;
80  fAlpha = alpha;
81  fTheta = theta;
82  fPhi = phi;
83  fTxy = TMath::Tan(alpha*TMath::DegToRad());
84  Double_t tth = TMath::Tan(theta*TMath::DegToRad());
85  Double_t ph = phi*TMath::DegToRad();
86  fTxz = tth*TMath::Cos(ph);
87  fTyz = tth*TMath::Sin(ph);
88  if ((fX<0) || (fY<0) || (fZ<0)) {
89 // printf("para : %f %f %f\n", fX, fY, fZ);
91  }
92  else ComputeBBox();
93 }
94 
95 ////////////////////////////////////////////////////////////////////////////////
96 /// Default constructor specifying minimum and maximum radius
97 
98 TGeoPara::TGeoPara(const char *name, Double_t dx, Double_t dy, Double_t dz, Double_t alpha,
99  Double_t theta, Double_t phi)
100  :TGeoBBox(name, 0, 0, 0)
101 {
103  fX = dx;
104  fY = dy;
105  fZ = dz;
106  fAlpha = alpha;
107  fTheta = theta;
108  fPhi = phi;
109  fTxy = TMath::Tan(alpha*TMath::DegToRad());
110  Double_t tth = TMath::Tan(theta*TMath::DegToRad());
111  Double_t ph = phi*TMath::DegToRad();
112  fTxz = tth*TMath::Cos(ph);
113  fTyz = tth*TMath::Sin(ph);
114  if ((fX<0) || (fY<0) || (fZ<0)) {
115 // printf("para : %f %f %f\n", fX, fY, fZ);
117  }
118  else ComputeBBox();
119 }
120 
121 ////////////////////////////////////////////////////////////////////////////////
122 /// Default constructor
123 /// - param[0] = dx
124 /// - param[1] = dy
125 /// - param[2] = dz
126 /// - param[3] = alpha
127 /// - param[4] = theta
128 /// - param[5] = phi
129 
131  :TGeoBBox(0, 0, 0)
132 {
134  SetDimensions(param);
135  if ((fX<0) || (fY<0) || (fZ<0)) SetShapeBit(kGeoRunTimeShape);
136  else ComputeBBox();
137 }
138 
139 ////////////////////////////////////////////////////////////////////////////////
140 /// destructor
141 
143 {
144 }
145 
146 ////////////////////////////////////////////////////////////////////////////////
147 /// Computes capacity of the shape in [length^3]
148 
150 {
151  Double_t capacity = 8.*fX*fY*fZ;
152  return capacity;
153 }
154 
155 ////////////////////////////////////////////////////////////////////////////////
156 /// compute bounding box
157 
159 {
161  Double_t dy = fY+fZ*TMath::Abs(fTyz);
162  Double_t dz = fZ;
163  TGeoBBox::SetBoxDimensions(dx, dy, dz);
164  memset(fOrigin, 0, 3*sizeof(Double_t));
165 }
166 
167 ////////////////////////////////////////////////////////////////////////////////
168 /// Compute normal to closest surface from POINT.
169 
170 void TGeoPara::ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
171 {
172  Double_t saf[3];
173  // distance from point to higher Z face
174  saf[0] = TMath::Abs(fZ-TMath::Abs(point[2])); // Z
175 
176  Double_t yt = point[1]-fTyz*point[2];
177  saf[1] = TMath::Abs(fY-TMath::Abs(yt)); // Y
178  // cos of angle YZ
179  Double_t cty = 1.0/TMath::Sqrt(1.0+fTyz*fTyz);
180 
181  Double_t xt = point[0]-fTxz*point[2]-fTxy*yt;
182  saf[2] = TMath::Abs(fX-TMath::Abs(xt)); // X
183  // cos of angle XZ
184  Double_t ctx = 1.0/TMath::Sqrt(1.0+fTxy*fTxy+fTxz*fTxz);
185  saf[2] *= ctx;
186  saf[1] *= cty;
187  Int_t i = TMath::LocMin(3,saf);
188  switch (i) {
189  case 0:
190  norm[0] = norm[1] = 0;
191  norm[2] = TMath::Sign(1.,dir[2]);
192  return;
193  case 1:
194  norm[0] = 0;
195  norm[1] = cty;
196  norm[2] = - fTyz*cty;
197  break;
198  case 2:
201  norm[2] = -TMath::Sin(fTheta*TMath::DegToRad());
202  }
203  if (norm[0]*dir[0]+norm[1]*dir[1]+norm[2]*dir[2]<0) {
204  norm[0] = -norm[0];
205  norm[1] = -norm[1];
206  norm[2] = -norm[2];
207  }
208 }
209 
210 ////////////////////////////////////////////////////////////////////////////////
211 /// test if point is inside this sphere
212 /// test Z range
213 
214 Bool_t TGeoPara::Contains(const Double_t *point) const
215 {
216  if (TMath::Abs(point[2]) > fZ) return kFALSE;
217  // check X and Y
218  Double_t yt=point[1]-fTyz*point[2];
219  if (TMath::Abs(yt) > fY) return kFALSE;
220  Double_t xt=point[0]-fTxz*point[2]-fTxy*yt;
221  if (TMath::Abs(xt) > fX) return kFALSE;
222  return kTRUE;
223 }
224 
225 ////////////////////////////////////////////////////////////////////////////////
226 /// compute distance from inside point to surface of the para
227 /// Boundary safe algorithm.
228 
229 Double_t TGeoPara::DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
230 {
231  if (iact<3 && safe) {
232  // compute safety
233  *safe = Safety(point, kTRUE);
234  if (iact==0) return TGeoShape::Big();
235  if (iact==1 && step<*safe) return TGeoShape::Big();
236  }
237  Double_t saf[2];
238  Double_t snxt = TGeoShape::Big();
239  Double_t s;
240  saf[0] = fZ+point[2];
241  saf[1] = fZ-point[2];
242  if (!TGeoShape::IsSameWithinTolerance(dir[2],0)) {
243  s = (dir[2]>0)?(saf[1]/dir[2]):(-saf[0]/dir[2]);
244  if (s<0) return 0.0;
245  if (s<snxt) snxt = s;
246  }
247  // distance from point to center axis on Y
248  Double_t yt = point[1]-fTyz*point[2];
249  saf[0] = fY+yt;
250  saf[1] = fY-yt;
251  Double_t dy = dir[1]-fTyz*dir[2];
253  s = (dy>0)?(saf[1]/dy):(-saf[0]/dy);
254  if (s<0) return 0.0;
255  if (s<snxt) snxt = s;
256  }
257  // distance from point to center axis on X
258  Double_t xt = point[0]-fTxz*point[2]-fTxy*yt;
259  saf[0] = fX+xt;
260  saf[1] = fX-xt;
261  Double_t dx = dir[0]-fTxz*dir[2]-fTxy*dy;
263  s = (dx>0)?(saf[1]/dx):(-saf[0]/dx);
264  if (s<0) return 0.0;
265  if (s<snxt) snxt = s;
266  }
267  return snxt;
268 }
269 
270 ////////////////////////////////////////////////////////////////////////////////
271 /// compute distance from inside point to surface of the para
272 
273 Double_t TGeoPara::DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact, Double_t step, Double_t *safe) const
274 {
275  if (iact<3 && safe) {
276  // compute safe distance
277  *safe = Safety(point, kFALSE);
278  if (iact==0) return TGeoShape::Big();
279  if (iact==1 && step<*safe) return TGeoShape::Big();
280  }
281  Bool_t in = kTRUE;
282  Double_t safz;
283  safz = TMath::Abs(point[2])-fZ;
284  if (safz>0) {
285  // outside Z
286  if (point[2]*dir[2]>=0) return TGeoShape::Big();
287  in = kFALSE;
288  }
289  Double_t yt=point[1]-fTyz*point[2];
290  Double_t safy = TMath::Abs(yt)-fY;
291  Double_t dy=dir[1]-fTyz*dir[2];
292  if (safy>0) {
293  if (yt*dy>=0) return TGeoShape::Big();
294  in = kFALSE;
295  }
296  Double_t xt=point[0]-fTxy*yt-fTxz*point[2];
297  Double_t safx = TMath::Abs(xt)-fX;
298  Double_t dx=dir[0]-fTxy*dy-fTxz*dir[2];
299  if (safx>0) {
300  if (xt*dx>=0) return TGeoShape::Big();
301  in = kFALSE;
302  }
303  // protection in case point is actually inside
304  if (in) {
305  if (safz>safx && safz>safy) {
306  if (point[2]*dir[2]>0) return TGeoShape::Big();
307  return 0.0;
308  }
309  if (safx>safy) {
310  if (xt*dx>0) return TGeoShape::Big();
311  return 0.0;
312  }
313  if (yt*dy>0) return TGeoShape::Big();
314  return 0.0;
315  }
316  Double_t xnew,ynew,znew;
317  if (safz>0) {
318  Double_t snxt = safz/TMath::Abs(dir[2]);
319  xnew = point[0]+snxt*dir[0];
320  ynew = point[1]+snxt*dir[1];
321  znew = (point[2]>0)?fZ:(-fZ);
322  Double_t ytn = ynew-fTyz*znew;
323  if (TMath::Abs(ytn)<=fY) {
324  Double_t xtn = xnew-fTxy*ytn-fTxz*znew;
325  if (TMath::Abs(xtn)<=fX) return snxt;
326  }
327  }
328  if (safy>0) {
329  Double_t snxt = safy/TMath::Abs(dy);
330  znew = point[2]+snxt*dir[2];
331  if (TMath::Abs(znew)<=fZ) {
332  Double_t ytn = (yt>0)?fY:(-fY);
333  xnew = point[0]+snxt*dir[0];
334  Double_t xtn = xnew-fTxy*ytn-fTxz*znew;
335  if (TMath::Abs(xtn)<=fX) return snxt;
336  }
337  }
338  if (safx>0) {
339  Double_t snxt = safx/TMath::Abs(dx);
340  znew = point[2]+snxt*dir[2];
341  if (TMath::Abs(znew)<=fZ) {
342  ynew = point[1]+snxt*dir[1];
343  Double_t ytn = ynew-fTyz*znew;
344  if (TMath::Abs(ytn)<=fY) return snxt;
345  }
346  }
347  return TGeoShape::Big();
348 }
349 
350 ////////////////////////////////////////////////////////////////////////////////
351 /// Divide this parallelepiped shape belonging to volume "voldiv" into ndiv equal volumes
352 /// called divname, from start position with the given step. Returns pointer
353 /// to created division cell volume. In case a wrong division axis is supplied,
354 /// returns pointer to volume to be divided.
355 
356 TGeoVolume *TGeoPara::Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv,
357  Double_t start, Double_t step)
358 {
359  TGeoShape *shape; //--- shape to be created
360  TGeoVolume *vol; //--- division volume to be created
361  TGeoVolumeMulti *vmulti; //--- generic divided volume
362  TGeoPatternFinder *finder; //--- finder to be attached
363  TString opt = ""; //--- option to be attached
364  Double_t end=start+ndiv*step;
365  switch (iaxis) {
366  case 1: //--- divide on X
367  shape = new TGeoPara(step/2, fY, fZ,fAlpha,fTheta, fPhi);
368  finder = new TGeoPatternParaX(voldiv, ndiv, start, end);
369  opt = "X";
370  break;
371  case 2: //--- divide on Y
372  shape = new TGeoPara(fX, step/2, fZ, fAlpha, fTheta, fPhi);
373  finder = new TGeoPatternParaY(voldiv, ndiv, start, end);
374  opt = "Y";
375  break;
376  case 3: //--- divide on Z
377  shape = new TGeoPara(fX, fY, step/2, fAlpha, fTheta, fPhi);
378  finder = new TGeoPatternParaZ(voldiv, ndiv, start, end);
379  opt = "Z";
380  break;
381  default:
382  Error("Divide", "Wrong axis type for division");
383  return 0;
384  }
385  vol = new TGeoVolume(divname, shape, voldiv->GetMedium());
386  vmulti = gGeoManager->MakeVolumeMulti(divname, voldiv->GetMedium());
387  vmulti->AddVolume(vol);
388  voldiv->SetFinder(finder);
389  finder->SetDivIndex(voldiv->GetNdaughters());
390  for (Int_t ic=0; ic<ndiv; ic++) {
391  voldiv->AddNodeOffset(vol, ic, start+step/2.+ic*step, opt.Data());
392  ((TGeoNodeOffset*)voldiv->GetNodes()->At(voldiv->GetNdaughters()-1))->SetFinder(finder);
393  }
394  return vmulti;
395 }
396 
397 ////////////////////////////////////////////////////////////////////////////////
398 /// Get range of shape for a given axis.
399 
401 {
402  xlo = 0;
403  xhi = 0;
404  Double_t dx = 0;
405  switch (iaxis) {
406  case 1:
407  xlo = -fX;
408  xhi = fX;
409  dx = xhi-xlo;
410  return dx;
411  case 2:
412  xlo = -fY;
413  xhi = fY;
414  dx = xhi-xlo;
415  return dx;
416  case 3:
417  xlo = -fZ;
418  xhi = fZ;
419  dx = xhi-xlo;
420  return dx;
421  }
422  return dx;
423 }
424 
425 ////////////////////////////////////////////////////////////////////////////////
426 /// Fill vector param[4] with the bounding cylinder parameters. The order
427 /// is the following : Rmin, Rmax, Phi1, Phi2
428 
430 {
432 }
433 
434 ////////////////////////////////////////////////////////////////////////////////
435 /// Fills real parameters of a positioned box inside this. Returns 0 if successful.
436 
437 Int_t TGeoPara::GetFittingBox(const TGeoBBox *parambox, TGeoMatrix *mat, Double_t &dx, Double_t &dy, Double_t &dz) const
438 {
439  dx=dy=dz=0;
440  if (mat->IsRotation()) {
441  Error("GetFittingBox", "cannot handle parametrized rotated volumes");
442  return 1; // ### rotation not accepted ###
443  }
444  //--> translate the origin of the parametrized box to the frame of this box.
445  Double_t origin[3];
446  mat->LocalToMaster(parambox->GetOrigin(), origin);
447  if (!Contains(origin)) {
448  Error("GetFittingBox", "wrong matrix - parametrized box is outside this");
449  return 1; // ### wrong matrix ###
450  }
451  //--> now we have to get the valid range for all parametrized axis
452  Double_t dd[3];
453  dd[0] = parambox->GetDX();
454  dd[1] = parambox->GetDY();
455  dd[2] = parambox->GetDZ();
456  //-> check if Z range is fixed
457  if (dd[2]<0) {
458  dd[2] = TMath::Min(origin[2]+fZ, fZ-origin[2]);
459  if (dd[2]<0) {
460  Error("GetFittingBox", "wrong matrix");
461  return 1;
462  }
463  }
464  if (dd[0]>=0 && dd[1]>=0) {
465  dx = dd[0];
466  dy = dd[1];
467  dz = dd[2];
468  return 0;
469  }
470  //-> check now range at Z = origin[2] +/- dd[2]
471  Double_t upper[8];
472  Double_t lower[8];
473  Double_t z=origin[2]-dd[2];
474  lower[0]=z*fTxz-fTxy*fY-fX;
475  lower[1]=-fY+z*fTyz;
476  lower[2]=z*fTxz+fTxy*fY-fX;
477  lower[3]=fY+z*fTyz;
478  lower[4]=z*fTxz+fTxy*fY+fX;
479  lower[5]=fY+z*fTyz;
480  lower[6]=z*fTxz-fTxy*fY+fX;
481  lower[7]=-fY+z*fTyz;
482  z=origin[2]+dd[2];
483  upper[0]=z*fTxz-fTxy*fY-fX;
484  upper[1]=-fY+z*fTyz;
485  upper[2]=z*fTxz+fTxy*fY-fX;
486  upper[3]=fY+z*fTyz;
487  upper[4]=z*fTxz+fTxy*fY+fX;
488  upper[5]=fY+z*fTyz;
489  upper[6]=z*fTxz-fTxy*fY+fX;
490  upper[7]=-fY+z*fTyz;
491 
492  for (Int_t iaxis=0; iaxis<2; iaxis++) {
493  if (dd[iaxis]>=0) continue;
494  Double_t ddmin = TGeoShape::Big();
495  for (Int_t ivert=0; ivert<4; ivert++) {
496  ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-lower[2*ivert+iaxis]));
497  ddmin = TMath::Min(ddmin, TMath::Abs(origin[iaxis]-upper[2*ivert+iaxis]));
498  }
499  dd[iaxis] = ddmin;
500  }
501  dx = dd[0];
502  dy = dd[1];
503  dz = dd[2];
504  return 0;
505 }
506 
507 ////////////////////////////////////////////////////////////////////////////////
508 /// in case shape has some negative parameters, these has to be computed
509 /// in order to fit the mother
510 
512 {
513  if (!TestShapeBit(kGeoRunTimeShape)) return 0;
514  if (!mother->TestShapeBit(kGeoPara)) {
515  Error("GetMakeRuntimeShape", "invalid mother");
516  return 0;
517  }
518  Double_t dx, dy, dz;
519  if (fX<0) dx=((TGeoPara*)mother)->GetX();
520  else dx=fX;
521  if (fY<0) dy=((TGeoPara*)mother)->GetY();
522  else dy=fY;
523  if (fZ<0) dz=((TGeoPara*)mother)->GetZ();
524  else dz=fZ;
525  return (new TGeoPara(dx, dy, dz, fAlpha, fTheta, fPhi));
526 }
527 
528 ////////////////////////////////////////////////////////////////////////////////
529 /// print shape parameters
530 
532 {
533  printf("*** Shape %s: TGeoPara ***\n", GetName());
534  printf(" dX = %11.5f\n", fX);
535  printf(" dY = %11.5f\n", fY);
536  printf(" dZ = %11.5f\n", fZ);
537  printf(" alpha = %11.5f\n", fAlpha);
538  printf(" theta = %11.5f\n", fTheta);
539  printf(" phi = %11.5f\n", fPhi);
540  printf(" Bounding box:\n");
542 }
543 
544 ////////////////////////////////////////////////////////////////////////////////
545 /// computes the closest distance from given point to this shape, according
546 /// to option. The matching point on the shape is stored in spoint.
547 
548 Double_t TGeoPara::Safety(const Double_t *point, Bool_t in) const
549 {
550  Double_t saf[3];
551  // distance from point to higher Z face
552  saf[0] = fZ-TMath::Abs(point[2]); // Z
553 
554  Double_t yt = point[1]-fTyz*point[2];
555  saf[1] = fY-TMath::Abs(yt); // Y
556  // cos of angle YZ
557  Double_t cty = 1.0/TMath::Sqrt(1.0+fTyz*fTyz);
558 
559  Double_t xt = point[0]-fTxz*point[2]-fTxy*yt;
560  saf[2] = fX-TMath::Abs(xt); // X
561  // cos of angle XZ
562  Double_t ctx = 1.0/TMath::Sqrt(1.0+fTxy*fTxy+fTxz*fTxz);
563  saf[2] *= ctx;
564  saf[1] *= cty;
565  if (in) return saf[TMath::LocMin(3,saf)];
566  for (Int_t i=0; i<3; i++) saf[i]=-saf[i];
567  return saf[TMath::LocMax(3,saf)];
568 }
569 
570 ////////////////////////////////////////////////////////////////////////////////
571 /// Save a primitive as a C++ statement(s) on output stream "out".
572 
573 void TGeoPara::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
574 {
575  if (TObject::TestBit(kGeoSavePrimitive)) return;
576  out << " // Shape: " << GetName() << " type: " << ClassName() << std::endl;
577  out << " dx = " << fX << ";" << std::endl;
578  out << " dy = " << fY << ";" << std::endl;
579  out << " dz = " << fZ << ";" << std::endl;
580  out << " alpha = " << fAlpha<< ";" << std::endl;
581  out << " theta = " << fTheta << ";" << std::endl;
582  out << " phi = " << fPhi << ";" << std::endl;
583  out << " TGeoShape *" << GetPointerName() << " = new TGeoPara(\"" << GetName() << "\",dx,dy,dz,alpha,theta,phi);" << std::endl;
585 }
586 
587 ////////////////////////////////////////////////////////////////////////////////
588 /// Set dimensions starting from an array.
589 
591 {
592  fX = param[0];
593  fY = param[1];
594  fZ = param[2];
595  fAlpha = param[3];
596  fTheta = param[4];
597  fPhi = param[5];
598  fTxy = TMath::Tan(param[3]*TMath::DegToRad());
599  Double_t tth = TMath::Tan(param[4]*TMath::DegToRad());
600  Double_t ph = param[5]*TMath::DegToRad();
601  fTxz = tth*TMath::Cos(ph);
602  fTyz = tth*TMath::Sin(ph);
603 }
604 
605 ////////////////////////////////////////////////////////////////////////////////
606 /// Create PARA mesh points
607 
609 {
610  if (!points) return;
611  Double_t txy = fTxy;
612  Double_t txz = fTxz;
613  Double_t tyz = fTyz;
614  *points++ = -fZ*txz-txy*fY-fX; *points++ = -fY-fZ*tyz; *points++ = -fZ;
615  *points++ = -fZ*txz+txy*fY-fX; *points++ = +fY-fZ*tyz; *points++ = -fZ;
616  *points++ = -fZ*txz+txy*fY+fX; *points++ = +fY-fZ*tyz; *points++ = -fZ;
617  *points++ = -fZ*txz-txy*fY+fX; *points++ = -fY-fZ*tyz; *points++ = -fZ;
618  *points++ = +fZ*txz-txy*fY-fX; *points++ = -fY+fZ*tyz; *points++ = +fZ;
619  *points++ = +fZ*txz+txy*fY-fX; *points++ = +fY+fZ*tyz; *points++ = +fZ;
620  *points++ = +fZ*txz+txy*fY+fX; *points++ = +fY+fZ*tyz; *points++ = +fZ;
621  *points++ = +fZ*txz-txy*fY+fX; *points++ = -fY+fZ*tyz; *points++ = +fZ;
622 }
623 
624 ////////////////////////////////////////////////////////////////////////////////
625 /// create sphere mesh points
626 
628 {
629  if (!points) return;
630  Double_t txy = fTxy;
631  Double_t txz = fTxz;
632  Double_t tyz = fTyz;
633  *points++ = -fZ*txz-txy*fY-fX; *points++ = -fY-fZ*tyz; *points++ = -fZ;
634  *points++ = -fZ*txz+txy*fY-fX; *points++ = +fY-fZ*tyz; *points++ = -fZ;
635  *points++ = -fZ*txz+txy*fY+fX; *points++ = +fY-fZ*tyz; *points++ = -fZ;
636  *points++ = -fZ*txz-txy*fY+fX; *points++ = -fY-fZ*tyz; *points++ = -fZ;
637  *points++ = +fZ*txz-txy*fY-fX; *points++ = -fY+fZ*tyz; *points++ = +fZ;
638  *points++ = +fZ*txz+txy*fY-fX; *points++ = +fY+fZ*tyz; *points++ = +fZ;
639  *points++ = +fZ*txz+txy*fY+fX; *points++ = +fY+fZ*tyz; *points++ = +fZ;
640  *points++ = +fZ*txz-txy*fY+fX; *points++ = -fY+fZ*tyz; *points++ = +fZ;
641 }
642 
643 ////////////////////////////////////////////////////////////////////////////////
644 /// fill size of this 3-D object
645 
646 void TGeoPara::Sizeof3D() const
647 {
649 }
650 
651 ////////////////////////////////////////////////////////////////////////////////
652 /// Check the inside status for each of the points in the array.
653 /// Input: Array of point coordinates + vector size
654 /// Output: Array of Booleans for the inside of each point
655 
656 void TGeoPara::Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const
657 {
658  for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
659 }
660 
661 ////////////////////////////////////////////////////////////////////////////////
662 /// Compute the normal for an array o points so that norm.dot.dir is positive
663 /// Input: Arrays of point coordinates and directions + vector size
664 /// Output: Array of normal directions
665 
666 void TGeoPara::ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize)
667 {
668  for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
669 }
670 
671 ////////////////////////////////////////////////////////////////////////////////
672 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
673 
674 void TGeoPara::DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
675 {
676  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
677 }
678 
679 ////////////////////////////////////////////////////////////////////////////////
680 /// Compute distance from array of input points having directions specified by dirs. Store output in dists
681 
682 void TGeoPara::DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step) const
683 {
684  for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
685 }
686 
687 ////////////////////////////////////////////////////////////////////////////////
688 /// Compute safe distance from each of the points in the input array.
689 /// Input: Array of point coordinates, array of statuses for these points, size of the arrays
690 /// Output: Safety values
691 
692 void TGeoPara::Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const
693 {
694  for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);
695 }
TGeoVolume::SetFinder
void SetFinder(TGeoPatternFinder *finder)
Definition: TGeoVolume.h:231
TGeoPara
Parallelepiped class.
Definition: TGeoPara.h:18
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
TGeoPara::GetAxisRange
virtual Double_t GetAxisRange(Int_t iaxis, Double_t &xlo, Double_t &xhi) const
Get range of shape for a given axis.
Definition: TGeoPara.cxx:400
TObject::TestBit
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
TGeoPara::~TGeoPara
virtual ~TGeoPara()
destructor
Definition: TGeoPara.cxx:142
TGeoVolume::GetNdaughters
Int_t GetNdaughters() const
Definition: TGeoVolume.h:349
Option_t
const char Option_t
Definition: RtypesCore.h:66
TGeoPara::ComputeBBox
virtual void ComputeBBox()
compute bounding box
Definition: TGeoPara.cxx:158
TGeoPara::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: TGeoPara.cxx:511
TMath::Cos
Double_t Cos(Double_t)
Definition: TMath.h:643
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
TGeoPara::Capacity
virtual Double_t Capacity() const
Computes capacity of the shape in [length^3].
Definition: TGeoPara.cxx:149
TMath::Sqrt
Double_t Sqrt(Double_t x)
Definition: TMath.h:691
TMath::DegToRad
constexpr Double_t DegToRad()
Conversion from degree to radian:
Definition: TMath.h:81
TMath::Tan
Double_t Tan(Double_t)
Definition: TMath.h:647
TGeoBBox::fOrigin
Double_t fOrigin[3]
Definition: TGeoBBox.h:24
Float_t
float Float_t
Definition: RtypesCore.h:57
TGeoPara::SavePrimitive
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TGeoPara.cxx:573
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:162
TGeoPatternParaY
Definition: TGeoPatternFinder.h:251
TGeoVolume.h
TGeoBBox::GetDZ
virtual Double_t GetDZ() const
Definition: TGeoBBox.h:76
TGeoPara::fZ
Double_t fZ
Definition: TGeoPara.h:23
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
TGeoPara::Contains
virtual Bool_t Contains(const Double_t *point) const
test if point is inside this sphere test Z range
Definition: TGeoPara.cxx:214
TMath::Abs
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
TGeoPara::Divide
virtual TGeoVolume * Divide(TGeoVolume *voldiv, const char *divname, Int_t iaxis, Int_t ndiv, Double_t start, Double_t step)
Divide this parallelepiped shape belonging to volume "voldiv" into ndiv equal volumes called divname,...
Definition: TGeoPara.cxx:356
TObjArray::At
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
TGeoPara::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: TGeoPara.cxx:682
TGeoPara::InspectShape
virtual void InspectShape() const
print shape parameters
Definition: TGeoPara.cxx:531
TGeoPara::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: TGeoPara.cxx:548
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
TGeoPara::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: TGeoPara.cxx:666
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
TGeoPara::Sizeof3D
virtual void Sizeof3D() const
fill size of this 3-D object
Definition: TGeoPara.cxx:646
TGeoVolumeMulti::AddVolume
void AddVolume(TGeoVolume *vol)
Add a volume with valid shape to the list of volumes.
Definition: TGeoVolume.cxx:2420
TGeoPara::fAlpha
Double_t fAlpha
Definition: TGeoPara.h:24
TGeoPatternParaZ
Definition: TGeoPatternFinder.h:287
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
TGeoShape::kGeoPara
@ kGeoPara
Definition: TGeoShape.h:45
TGeoShape
Base abstract class for all shapes.
Definition: TGeoShape.h:26
TGeoBBox::InspectShape
virtual void InspectShape() const
Prints shape parameters.
Definition: TGeoBBox.cxx:790
TObject::SetBit
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
TGeoPara::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 inside point to surface of the para
Definition: TGeoPara.cxx:273
TGeoPara::fY
Double_t fY
Definition: TGeoPara.h:22
TGeoVolume::GetNodes
TObjArray * GetNodes()
Definition: TGeoVolume.h:167
TGeoBBox::GetBoundingCylinder
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoBBox.cxx:580
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
TGeoPara::SetPoints
virtual void SetPoints(Double_t *points) const
Create PARA mesh points.
Definition: TGeoPara.cxx:608
TGeoShape::kGeoSavePrimitive
@ kGeoSavePrimitive
Definition: TGeoShape.h:65
TGeoPara::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: TGeoPara.cxx:656
TGeoPara::SetDimensions
virtual void SetDimensions(Double_t *param)
Set dimensions starting from an array.
Definition: TGeoPara.cxx:590
TGeoBBox::SetBoxDimensions
void SetBoxDimensions(Double_t dx, Double_t dy, Double_t dz, Double_t *origin=0)
Set parameters of the box.
Definition: TGeoBBox.cxx:897
TGeoPara::fX
Double_t fX
Definition: TGeoPara.h:21
TMath::Sin
Double_t Sin(Double_t)
Definition: TMath.h:639
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
TGeoBBox::GetOrigin
virtual const Double_t * GetOrigin() const
Definition: TGeoBBox.h:77
TGeoBBox
Box class.
Definition: TGeoBBox.h:18
TGeoPara::fTyz
Double_t fTyz
Definition: TGeoPara.h:29
TGeoPara::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 para Boundary safe algorithm.
Definition: TGeoPara.cxx:229
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
TMath::Min
Short_t Min(Short_t a, Short_t b)
Definition: TMathBase.h:180
TGeoMatrix.h
TGeoManager.h
TGeoPara::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: TGeoPara.cxx:437
TGeoPara::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: TGeoPara.cxx:692
TGeoPara::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: TGeoPara.cxx:674
TGeoShape::GetPointerName
const char * GetPointerName() const
Provide a pointer name containing uid.
Definition: TGeoShape.cxx:699
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
TGeoBBox::GetDY
virtual Double_t GetDY() const
Definition: TGeoBBox.h:75
TGeoVolumeMulti
Volume families.
Definition: TGeoVolume.h:254
TGeoPara::ComputeNormal
virtual void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm)
Compute normal to closest surface from POINT.
Definition: TGeoPara.cxx:170
TGeoPara::GetBoundingCylinder
virtual void GetBoundingCylinder(Double_t *param) const
Fill vector param[4] with the bounding cylinder parameters.
Definition: TGeoPara.cxx:429
points
point * points
Definition: X3DBuffer.c:22
TGeoShape::IsSameWithinTolerance
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
Definition: TGeoShape.cxx:326
TGeoPara::fPhi
Double_t fPhi
Definition: TGeoPara.h:26
name
char name[80]
Definition: TGX11.cxx:110
TGeoPara.h
TGeoPara::TGeoPara
TGeoPara()
Default constructor.
Definition: TGeoPara.cxx:57
TGeoPatternParaX
Definition: TGeoPatternFinder.h:218
TGeoPara::fTxz
Double_t fTxz
Definition: TGeoPara.h:28
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
TGeoPara::fTheta
Double_t fTheta
Definition: TGeoPara.h:25
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
TGeoPara::fTxy
Double_t fTxy
Definition: TGeoPara.h:27
TGeoVolume
TGeoVolume, TGeoVolumeMulti, TGeoVolumeAssembly are the volume classes.
Definition: TGeoVolume.h:49
TMath.h
int