Geant4  10.02.p02
G4VTwistSurface.icc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4VTwistSurface.icc 66356 2012-12-18 09:02:32Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------------
31 // G4VTwistSurface class inline methods
32 //
33 // Author:
34 // 01-Aug-2002 - Kotoyo Hoshina (hoshina@hepburn.s.chiba-u.ac.jp)
35 //
36 // History:
37 // 13-Nov-2003 - O.Link (Oliver.Link@cern.ch), Integration in Geant4
38 // from original version in Jupiter-2.5.02 application.
39 // --------------------------------------------------------------------
40 
41 //=====================================================================
42 //* DistanceToPlaneWithV ----------------------------------------------
43 
44 inline
45 G4double G4VTwistSurface::DistanceToPlaneWithV(const G4ThreeVector &p,
46  const G4ThreeVector &v,
47  const G4ThreeVector &x0,
48  const G4ThreeVector &n0,
49  G4ThreeVector &xx)
50 {
51  G4double q = n0 * v;
52  G4double t = kInfinity;
53  if (q) { t = (n0 * (x0 - p)) / q; }
54  xx = p + t * v;
55  return t;
56 }
57 
58 //=====================================================================
59 //* DistanceToPlane ---------------------------------------------------
60 
61 inline
62 G4double G4VTwistSurface::DistanceToPlane(const G4ThreeVector &p,
63  const G4ThreeVector &x0,
64  const G4ThreeVector &n0,
65  G4ThreeVector &xx)
66 {
67  // DistanceToPlane :
68  // Calculate distance to plane in local coordinate,
69  // then return distance and global intersection points.
70  //
71  // p - location of flying particle
72  // x0 - reference point of surface
73  // xx - a foot of perpendicular line from p to the plane
74  // t - distance from xx to p
75  // n - a unit normal of this plane from plane to p.
76  //
77  // equation of plane:
78  // n*(x - x0) = 0;
79  //
80  // vector to xx:
81  // xx = p - t*n
82  //
83  // where
84  // t = n * (p - x0) / std::fabs(n)
85  //
86  G4double t;
87  G4ThreeVector n = n0.unit();
88  t = n * (p - x0);
89  xx = p - t * n;
90  return t;
91 }
92 
93 //=====================================================================
94 //* DistanceToPlane ---------------------------------------------------
95 
96 inline
97 G4double G4VTwistSurface::DistanceToPlane(const G4ThreeVector &p,
98  const G4ThreeVector &x0,
99  const G4ThreeVector &t1,
100  const G4ThreeVector &t2,
101  G4ThreeVector &xx,
102  G4ThreeVector &n)
103 {
104  // DistanceToPlane :
105  // Calculate distance to plane in local coordinate,
106  // then return distance and global intersection points.
107  // t1 - 1st. vector lying on the plane
108  // t2 - 2nd. vector lying on the plane
109 
110  n = (t1.cross(t2)).unit();
111  return DistanceToPlane(p, x0, n, xx);
112 }
113 
114 //=====================================================================
115 //* DistanceToLine ----------------------------------------------------
116 
117 inline
118 G4double G4VTwistSurface::DistanceToLine(const G4ThreeVector &p,
119  const G4ThreeVector &x0,
120  const G4ThreeVector &d,
121  G4ThreeVector &xx)
122 {
123  // DistanceToLine :
124  // Calculate distance to line,
125  // then return distance and global intersection points.
126  //
127  // p - location of flying particle
128  // x0 - reference point of line
129  // d - direction vector of line
130  // xx - a foot of perpendicular line from p to the plane
131  // t - distance from xx to p
132  //
133  // Equation
134  //
135  // distance^2 = |(xx - p)|^2
136  // with
137  // xx = x0 + t*d
138  //
139  // (d/dt)distance^2 = (d/dt)|((x0 + t*d) - p)|^2
140  // = 2*t*|d|^2 + 2*d*(x0 - p)
141  // = 0 // smallest distance
142  // then
143  // t = - d*(x0 - p) / |d|^2
144  //
145 
146  G4double t;
147  G4ThreeVector dir = d.unit();
148  t = - dir * (x0 - p); // |dir|^2 = 1.
149  xx = x0 + t * dir;
150 
151  G4ThreeVector dist = xx - p;
152  return dist.mag();
153 }
154 
155 //=====================================================================
156 //* IsAxis0 -----------------------------------------------------------
157 
158 inline
159 G4bool G4VTwistSurface::IsAxis0(G4int areacode) const
160 {
161  if (areacode & sAxis0) return true;
162  return false;
163 }
164 
165 //=====================================================================
166 //* IsAxis1 -----------------------------------------------------------
167 
168 inline
169 G4bool G4VTwistSurface::IsAxis1(G4int areacode) const
170 {
171  if (areacode & sAxis1) return true;
172  return false;
173 }
174 
175 //=====================================================================
176 //* IsOutside ---------------------------------------------------------
177 
178 inline
179 G4bool G4VTwistSurface::IsOutside(G4int areacode) const
180 {
181  if (areacode & sInside) return false;
182  return true;
183 }
184 
185 //=====================================================================
186 //* IsInside ----------------------------------------------------------
187 
188 inline
189 G4bool G4VTwistSurface::IsInside(G4int areacode, G4bool testbitmode) const
190 {
191  if (areacode & sInside) {
192  if (testbitmode) {
193  return true;
194  } else {
195  if (!((areacode & sBoundary) || (areacode & sCorner))) return true;
196  }
197  }
198  return false;
199 }
200 
201 //=====================================================================
202 //* IsBoundary --------------------------------------------------------
203 
204 inline
205 G4bool G4VTwistSurface::IsBoundary(G4int areacode, G4bool testbitmode) const
206 {
207  if ((areacode & sBoundary) == sBoundary) {
208  if (testbitmode) {
209  return true;
210  } else {
211  if ((areacode & sInside) == sInside) return true;
212  }
213  }
214  return false;
215 }
216 
217 //=====================================================================
218 //* IsCorner ----------------------------------------------------------
219 
220 inline
221 G4bool G4VTwistSurface::IsCorner(G4int areacode, G4bool testbitmode) const
222 {
223  if ((areacode & sCorner) == sCorner) {
224  if (testbitmode) {
225  return true;
226  } else {
227  if ((areacode & sInside) == sInside) return true;
228  }
229  }
230  return false;
231 }
232 
233 //=====================================================================
234 //* GetAxisType -------------------------------------------------------
235 
236 inline
237 G4int G4VTwistSurface::GetAxisType(G4int areacode, G4int whichaxis) const
238 {
239  G4int axiscode = areacode & sAxisMask & whichaxis;
240 
241  if (axiscode == (sAxisX & sAxis0) ||
242  axiscode == (sAxisX & sAxis1)) {
243  return sAxisX;
244  } else if (axiscode == (sAxisY & sAxis0) ||
245  axiscode == (sAxisY & sAxis1)) {
246  return sAxisY;
247  } else if (axiscode == (sAxisZ & sAxis0) ||
248  axiscode == (sAxisZ & sAxis1)) {
249  return sAxisZ;
250  } else if (axiscode == (sAxisRho & sAxis0) ||
251  axiscode == (sAxisRho & sAxis1)) {
252  return sAxisRho;
253  } else if (axiscode == (sAxisPhi & sAxis0) ||
254  axiscode == (sAxisPhi & sAxis1)) {
255  return sAxisPhi;
256  } else {
257  std::ostringstream message;
258  message << "Configuration not supported." << G4endl
259  << " areacode = " << areacode;
260  G4Exception("G4VTwistSurface::GetAxisType()","GeomSolids0001",
261  FatalException, message);
262  }
263  return 1;
264 }
265 
266 //=====================================================================
267 //* ComputeGlobalPoint ------------------------------------------------
268 
269 inline
270 G4ThreeVector G4VTwistSurface::ComputeGlobalPoint(const G4ThreeVector &lp) const
271 {
272  return fRot * G4ThreeVector(lp) + fTrans;
273 }
274 
275 //=====================================================================
276 //* ComputeGlobalPoint ------------------------------------------------
277 
278 inline
279 G4ThreeVector G4VTwistSurface::ComputeLocalPoint(const G4ThreeVector &gp) const
280 {
281  return fRot.inverse() * ( G4ThreeVector(gp) - fTrans ) ;
282 }
283 
284 //=====================================================================
285 //* ComputeGlobalDirection --------------------------------------------
286 
287 inline
288 G4ThreeVector G4VTwistSurface::ComputeGlobalDirection(const G4ThreeVector &lp) const
289 {
290  return fRot * G4ThreeVector(lp);
291 }
292 
293 //=====================================================================
294 //* ComputeLocalDirection ---------------------------------------------
295 
296 inline
297 G4ThreeVector G4VTwistSurface::ComputeLocalDirection(const G4ThreeVector &gp) const
298 {
299  return fRot.inverse() * G4ThreeVector(gp);
300 }
301 
302 //=====================================================================
303 //* SetNeighbours -----------------------------------------------------
304 
305 inline
306 void G4VTwistSurface::SetNeighbours(G4VTwistSurface* axis0min, G4VTwistSurface* axis1min,
307  G4VTwistSurface* axis0max, G4VTwistSurface* axis1max)
308 {
309  fNeighbours[0] = axis0min;
310  fNeighbours[1] = axis1min;
311  fNeighbours[2] = axis0max;
312  fNeighbours[3] = axis1max;
313 }
314 
315 //=====================================================================
316 //* GetNeighbours -----------------------------------------------------
317 
318 inline
319 G4int G4VTwistSurface::GetNeighbours(G4int areacode, G4VTwistSurface** surfaces)
320 {
321 
322  int sAxis0Min = sAxis0 & sAxisMin ;
323  int sAxis1Min = sAxis1 & sAxisMin ;
324  int sAxis0Max = sAxis0 & sAxisMax ;
325  int sAxis1Max = sAxis1 & sAxisMax ;
326 
327  G4int i = 0;
328 
329  if ( (areacode & sAxis0Min ) == sAxis0Min ) {
330  surfaces[i] = fNeighbours[0] ;
331  i++ ;
332  }
333 
334  if ( ( areacode & sAxis1Min ) == sAxis1Min ) {
335  surfaces[i] = fNeighbours[1] ;
336  i++ ;
337  if ( i == 2 ) return i ;
338  }
339 
340  if ( ( areacode & sAxis0Max ) == sAxis0Max ) {
341  surfaces[i] = fNeighbours[2] ;
342  i++ ;
343  if ( i == 2 ) return i ;
344  }
345 
346  if ( ( areacode & sAxis1Max ) == sAxis1Max ) {
347  surfaces[i] = fNeighbours[3] ;
348  i++ ;
349  if ( i == 2 ) return i ;
350  }
351 
352  return i ;
353 }
354 
355 //=====================================================================
356 //* GetCorner ---------------------------------------------------------
357 
358 inline
359 G4ThreeVector G4VTwistSurface::GetCorner(G4int areacode) const
360 {
361  if (!(areacode & sCorner)){
362  std::ostringstream message;
363  message << "Area code must represent corner." << G4endl
364  << " areacode = " << areacode;
365  G4Exception("G4VTwistSurface::GetCorner()","GeomSolids0002",
366  FatalException, message);
367  }
368 
369  if ((areacode & sC0Min1Min) == sC0Min1Min) {
370  return fCorners[0];
371  } else if ((areacode & sC0Max1Min) == sC0Max1Min) {
372  return fCorners[1];
373  } else if ((areacode & sC0Max1Max) == sC0Max1Max) {
374  return fCorners[2];
375  } else if ((areacode & sC0Min1Max) == sC0Min1Max) {
376  return fCorners[3];
377  } else {
378  std::ostringstream message;
379  message << "Configuration not supported." << G4endl
380  << " areacode = " << areacode;
381  G4Exception("G4VTwistSurface::GetCorner()", "GeomSolids0001",
382  FatalException, message);
383  }
384  return fCorners[0];
385 }