EGSnrc C++ class library  Report PIRS-898 (2021)
Iwan Kawrakow, Ernesto Mainegra-Hing, Frederic Tessier, Reid Townson and Blake Walters
egs_circle_perpendicular.h
Go to the documentation of this file.
1 /*
2 ###############################################################################
3 #
4 # EGSnrc egs++ circle perpendicular shape headers
5 # Copyright (C) 2015 National Research Council Canada
6 #
7 # This file is part of EGSnrc.
8 #
9 # EGSnrc is free software: you can redistribute it and/or modify it under
10 # the terms of the GNU Affero General Public License as published by the
11 # Free Software Foundation, either version 3 of the License, or (at your
12 # option) any later version.
13 #
14 # EGSnrc is distributed in the hope that it will be useful, but WITHOUT ANY
15 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16 # FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for
17 # more details.
18 #
19 # You should have received a copy of the GNU Affero General Public License
20 # along with EGSnrc. If not, see <http://www.gnu.org/licenses/>.
21 #
22 ###############################################################################
23 #
24 # Author: Reid Townson, 2020
25 #
26 # Contributors:
27 #
28 ###############################################################################
29 */
30 
31 
37 #ifndef EGS_CIRCLE_PERPENDICULAR_
38 #define EGS_CIRCLE_PERPENDICULAR_
39 
40 #include "egs_shapes.h"
41 #include "egs_rndm.h"
42 #include "egs_math.h"
43 
44 #ifdef WIN32
45 
46  #ifdef BUILD_CIRCLE_PERPENDICULAR_DLL
47  #define EGS_CIRCLE_PERPENDICULAR_EXPORT __declspec(dllexport)
48  #else
49  #define EGS_CIRCLE_PERPENDICULAR_EXPORT __declspec(dllimport)
50  #endif
51  #define EGS_CIRCLE_PERPENDICULAR_LOCAL
52 
53 #else
54 
55  #ifdef HAVE_VISIBILITY
56  #define EGS_CIRCLE_PERPENDICULAR_EXPORT __attribute__ ((visibility ("default")))
57  #define EGS_CIRCLE_PERPENDICULAR_LOCAL __attribute__ ((visibility ("hidden")))
58  #else
59  #define EGS_CIRCLE_PERPENDICULAR_EXPORT
60  #define EGS_CIRCLE_PERPENDICULAR_LOCAL
61  #endif
62 
63 #endif
64 
109 class EGS_CIRCLE_PERPENDICULAR_EXPORT EGS_CirclePerpendicularShape : public EGS_SurfaceShape {
110 
111 public:
112 
115  EGS_CirclePerpendicularShape(EGS_Float Xo, EGS_Float Yo, EGS_Float R, EGS_Float R_i = 0,
116  const string &Name="",EGS_ObjectFactory *f=0) :
117  EGS_SurfaceShape(Name,f), xo(Xo), yo(Yo), ro(R_i), dr(R-R_i) {
118  otype = "circle";
119  if (dr < 0) {
120  ro = R;
121  dr = R_i - R;
122  }
123  A = M_PI*dr*(dr + 2*ro);
124  };
127  EGS_Float r = ro + dr*sqrt(rndm->getUniform());
128  EGS_Float cphi, sphi;
129  rndm->getAzimuth(cphi,sphi);
130  return EGS_Vector(xo + r*cphi, yo + r*sphi, 0);
131  };
132 
133  void getPointSourceDirection(const EGS_Vector &Xo,
134  EGS_RandomGenerator *rndm, EGS_Vector &u, EGS_Float &wt) {
135 
136  // Perform user requested transformations to a point on the source surface
137  // This is effectively the same as transforming the target position
138  // because we're just calculating the direction between the two
139  EGS_Vector xo = T ? Xo*(*T) : Xo;
140 
141  // Get a point on the circle target surface
142  EGS_Vector x = getPoint(rndm);
143  u = x - xo;
144  EGS_Float d2i = 1/u.length2(), di = sqrt(d2i);
145  u *= di;
146 
147  // Calculate the angle between the normal to the circle surface and the u vector between points
148  EGS_Vector perpToCircle = EGS_Vector(0, 0, 1);
149  EGS_Float angleBetween = std::acos(u * perpToCircle / (perpToCircle.length() * u.length()));
150 
151  // We will rotate about a vector perpendicular to u and the target surface normal
152  // Check against fabs(u.z) to account for both parallel and anti-parallel cases
153  EGS_Vector rotateAbout;
154  if ((fabs(u.z) - perpToCircle.z) < epsilon) {
155  rotateAbout = perpToCircle;
156  }
157  else {
158  rotateAbout = u.times(perpToCircle);
159  }
160  EGS_RotationMatrix rotation = EGS_RotationMatrix::rotV(-angleBetween, rotateAbout);
161  EGS_AffineTransform *transform = new EGS_AffineTransform(rotation);
162  if (transform) {
163  // Transform the point on the target surface by this rotation
164  transform->rotate(x);
165 
166  // Calculate the new u vector
167  u = x - xo;
168  }
169  delete transform;
170 
171  d2i = 1/u.length2(), di = sqrt(d2i);
172  u *= di;
173  wt = A*u.length()*d2i;
174  if (T) {
175  T->rotate(u);
176  }
177  };
178 
179 protected:
180 
181  EGS_Float xo, yo, ro, dr;
182 };
183 
184 #endif
A surface shape.
Definition: egs_shapes.h:252
A circle shape perpendicular to source particles.
A class providing affine transformations.
A class for vector rotations.
A class representing 3D vectors.
Definition: egs_vector.h:56
void getAzimuth(EGS_Float &cphi, EGS_Float &sphi)
Sets cphi and sphi to the cosine and sine of a random angle uniformely distributed between 0 and ...
Definition: egs_rndm.h:138
EGS_CirclePerpendicularShape(EGS_Float Xo, EGS_Float Yo, EGS_Float R, EGS_Float R_i=0, const string &Name="", EGS_ObjectFactory *f=0)
Conctruct a circle with midpoint given by Xo and Yo, radius R and innder radius R_i.
virtual EGS_Vector getPoint(EGS_RandomGenerator *rndm)
Sample and return a random 3D vector.
Definition: egs_shapes.h:149
Base random number generator class. All random number generators should be derived from this class...
Definition: egs_rndm.h:67
EGS_Float z
z-component
Definition: egs_vector.h:62
EGS_RandomGenerator class header file.
void rotate(EGS_Vector &v) const
Applies the rotation to the vector v.
EGS_BaseShape and shape classes header file.
An object factory.
static EGS_RotationMatrix rotV(EGS_Float phi, const EGS_Vector &v)
Returns a rotation by the angle phi around the axis defined by the vector v.
Attempts to fix broken math header files.
EGS_Float getUniform()
Returns a random number uniformly distributed between zero (inclusive) and 1 (exclusive).
Definition: egs_rndm.h:103
const EGS_Float epsilon
The epsilon constant for floating point comparisons.
Definition: egs_functions.h:61
void getPointSourceDirection(const EGS_Vector &Xo, EGS_RandomGenerator *rndm, EGS_Vector &u, EGS_Float &wt)
Get a random direction given a source position Xo.
Definition: egs_shapes.h:275