| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration and of QinetiQ Ltd, *
20 // * subject to DEFCON 705 IPR conditions. 20 // * subject to DEFCON 705 IPR conditions. *
21 // * By using, copying, modifying or distri 21 // * By using, copying, modifying or distributing the software (or *
22 // * any work based on the software) you ag 22 // * any work based on the software) you agree to acknowledge its *
23 // * use in resulting scientific publicati 23 // * use in resulting scientific publications, and indicate your *
24 // * acceptance of all terms of the Geant4 Sof 24 // * acceptance of all terms of the Geant4 Software license. *
25 // ******************************************* 25 // ********************************************************************
26 // 26 //
27 // G4TessellatedGeometryAlgorithms implementat << 27 //
>> 28 // $Id: G4TessellatedGeometryAlgorithms.cc 66356 2012-12-18 09:02:32Z gcosmo $
>> 29 //
>> 30 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>> 31 //
>> 32 // CHANGE HISTORY
>> 33 // --------------
28 // 34 //
29 // 07 August 2007, P R Truscott, QinetiQ Ltd, 35 // 07 August 2007, P R Truscott, QinetiQ Ltd, UK - Created, with member
30 // functions based on the work 36 // functions based on the work of Rickard Holmberg.
>> 37 //
>> 38 // 26 September 2007
>> 39 // P R Truscott, qinetiQ Ltd, UK
>> 40 // Updated to assign values of location array, not update
>> 41 // just the pointer.
>> 42 //
31 // 12 October 2012, M Gayer, CERN, - Reviewed 43 // 12 October 2012, M Gayer, CERN, - Reviewed optimized implementation.
32 // ------------------------------------------- << 44 //
>> 45 ///////////////////////////////////////////////////////////////////////////////
33 46
34 #include "G4TessellatedGeometryAlgorithms.hh" 47 #include "G4TessellatedGeometryAlgorithms.hh"
35 48
36 #include <cfloat> 49 #include <cfloat>
37 50
38 ////////////////////////////////////////////// 51 ///////////////////////////////////////////////////////////////////////////////
39 // 52 //
40 // IntersectLineAndTriangle2D 53 // IntersectLineAndTriangle2D
41 // 54 //
42 // Determines whether there is an intersection 55 // Determines whether there is an intersection between a line defined
43 // by r = p + s.v and a triangle defined by ve 56 // by r = p + s.v and a triangle defined by verticies p0, p0+e0 and p0+e1.
44 // 57 //
45 // Here: 58 // Here:
46 // p = 2D vector 59 // p = 2D vector
47 // s = scaler on [0,infinity) 60 // s = scaler on [0,infinity)
48 // v = 2D vector 61 // v = 2D vector
49 // p0, e0 and e1 are 2D vectors 62 // p0, e0 and e1 are 2D vectors
50 // Information about where the intersection oc 63 // Information about where the intersection occurs is returned in the
51 // variable location. 64 // variable location.
52 // 65 //
53 // This is based on the work of Rickard Holmbe 66 // This is based on the work of Rickard Holmberg.
54 // 67 //
55 G4bool G4TessellatedGeometryAlgorithms::Inters 68 G4bool G4TessellatedGeometryAlgorithms::IntersectLineAndTriangle2D (
56 const G4TwoVector& p, const G4TwoVector& v, << 69 const G4TwoVector &p, const G4TwoVector &v,
57 const G4TwoVector& p0, const G4TwoVector& e0 << 70 const G4TwoVector &p0, const G4TwoVector &e0, const G4TwoVector &e1,
58 G4TwoVector location[2]) 71 G4TwoVector location[2])
59 { 72 {
60 G4TwoVector loc0[2]; 73 G4TwoVector loc0[2];
61 G4int e0i = IntersectLineAndLineSegment2D (p 74 G4int e0i = IntersectLineAndLineSegment2D (p,v,p0,e0,loc0);
62 if (e0i == 2) 75 if (e0i == 2)
63 { 76 {
64 location[0] = loc0[0]; 77 location[0] = loc0[0];
65 location[1] = loc0[1]; 78 location[1] = loc0[1];
66 return true; 79 return true;
67 } 80 }
68 81
69 G4TwoVector loc1[2]; 82 G4TwoVector loc1[2];
70 G4int e1i = IntersectLineAndLineSegment2D (p 83 G4int e1i = IntersectLineAndLineSegment2D (p,v,p0,e1,loc1);
71 if (e1i == 2) 84 if (e1i == 2)
72 { 85 {
73 location[0] = loc1[0]; 86 location[0] = loc1[0];
74 location[1] = loc1[1]; 87 location[1] = loc1[1];
75 return true; 88 return true;
76 } 89 }
77 90
78 if ((e0i == 1) && (e1i == 1)) 91 if ((e0i == 1) && (e1i == 1))
79 { 92 {
80 if ((loc0[0]-p).mag2() < (loc1[0]-p).mag2( 93 if ((loc0[0]-p).mag2() < (loc1[0]-p).mag2())
81 { 94 {
82 location[0] = loc0[0]; 95 location[0] = loc0[0];
83 location[1] = loc1[0]; 96 location[1] = loc1[0];
84 } 97 }
85 else 98 else
86 { 99 {
87 location[0] = loc1[0]; 100 location[0] = loc1[0];
88 location[1] = loc0[0]; 101 location[1] = loc0[0];
89 } 102 }
90 return true; 103 return true;
91 } 104 }
92 105
93 G4TwoVector p1 = p0 + e0; 106 G4TwoVector p1 = p0 + e0;
94 G4TwoVector DE = e1 - e0; 107 G4TwoVector DE = e1 - e0;
95 G4TwoVector loc2[2]; 108 G4TwoVector loc2[2];
96 G4int e2i = IntersectLineAndLineSegment2D (p 109 G4int e2i = IntersectLineAndLineSegment2D (p,v,p1,DE,loc2);
97 if (e2i == 2) 110 if (e2i == 2)
98 { 111 {
99 location[0] = loc2[0]; 112 location[0] = loc2[0];
100 location[1] = loc2[1]; 113 location[1] = loc2[1];
101 return true; 114 return true;
102 } 115 }
103 116
104 if ((e0i == 0) && (e1i == 0) && (e2i == 0)) 117 if ((e0i == 0) && (e1i == 0) && (e2i == 0)) return false;
105 118
106 if ((e0i == 1) && (e2i == 1)) 119 if ((e0i == 1) && (e2i == 1))
107 { 120 {
108 if ((loc0[0]-p).mag2() < (loc2[0]-p).mag2( 121 if ((loc0[0]-p).mag2() < (loc2[0]-p).mag2())
109 { 122 {
110 location[0] = loc0[0]; 123 location[0] = loc0[0];
111 location[1] = loc2[0]; 124 location[1] = loc2[0];
112 } 125 }
113 else 126 else
114 { 127 {
115 location[0] = loc2[0]; 128 location[0] = loc2[0];
116 location[1] = loc0[0]; 129 location[1] = loc0[0];
117 } 130 }
118 return true; 131 return true;
119 } 132 }
120 133
121 if ((e1i == 1) && (e2i == 1)) 134 if ((e1i == 1) && (e2i == 1))
122 { 135 {
123 if ((loc1[0]-p).mag2() < (loc2[0]-p).mag2( 136 if ((loc1[0]-p).mag2() < (loc2[0]-p).mag2())
124 { 137 {
125 location[0] = loc1[0]; 138 location[0] = loc1[0];
126 location[1] = loc2[0]; 139 location[1] = loc2[0];
127 } 140 }
128 else 141 else
129 { 142 {
130 location[0] = loc2[0]; 143 location[0] = loc2[0];
131 location[1] = loc1[0]; 144 location[1] = loc1[0];
132 } 145 }
133 return true; 146 return true;
134 } 147 }
135 148
136 return false; 149 return false;
137 } 150 }
138 151
139 ////////////////////////////////////////////// 152 ///////////////////////////////////////////////////////////////////////////////
140 // 153 //
141 // IntersectLineAndLineSegment2D 154 // IntersectLineAndLineSegment2D
142 // 155 //
143 // Determines whether there is an intersection 156 // Determines whether there is an intersection between a line defined
144 // by r = p0 + s.d0 and a line-segment with en 157 // by r = p0 + s.d0 and a line-segment with endpoints p1 and p1+d1.
145 // Here: 158 // Here:
146 // p0 = 2D vector 159 // p0 = 2D vector
147 // s = scaler on [0,infinity) 160 // s = scaler on [0,infinity)
148 // d0 = 2D vector 161 // d0 = 2D vector
149 // p1 and d1 are 2D vectors 162 // p1 and d1 are 2D vectors
150 // 163 //
151 // This function returns: 164 // This function returns:
152 // 0 - if there is no intersection; 165 // 0 - if there is no intersection;
153 // 1 - if there is a unique intersection; 166 // 1 - if there is a unique intersection;
154 // 2 - if the line and line-segments overlap, 167 // 2 - if the line and line-segments overlap, and the intersection is a
155 // segment itself. 168 // segment itself.
156 // Information about where the intersection oc 169 // Information about where the intersection occurs is returned in the
157 // as ??. 170 // as ??.
158 // 171 //
159 // This is based on the work of Rickard Holmbe 172 // This is based on the work of Rickard Holmberg as well as material published
160 // by Philip J Schneider and David H Eberly, " 173 // by Philip J Schneider and David H Eberly, "Geometric Tools for Computer
161 // Graphics," ISBN 1-55860-694-0, pp 244-245, 174 // Graphics," ISBN 1-55860-694-0, pp 244-245, 2003.
162 // 175 //
163 G4int G4TessellatedGeometryAlgorithms::Interse 176 G4int G4TessellatedGeometryAlgorithms::IntersectLineAndLineSegment2D (
164 const G4TwoVector& p0, const G4TwoVector& d0 << 177 const G4TwoVector &p0, const G4TwoVector &d0,
165 const G4TwoVector& p1, const G4TwoVector& d1 << 178 const G4TwoVector &p1, const G4TwoVector &d1,
>> 179 G4TwoVector location[2])
166 { 180 {
167 G4TwoVector e = p1 - p0; 181 G4TwoVector e = p1 - p0;
168 G4double kross = cross(d0,d1); 182 G4double kross = cross(d0,d1);
169 G4double sqrKross = kross * kross; 183 G4double sqrKross = kross * kross;
170 G4double sqrLen0 = d0.mag2(); 184 G4double sqrLen0 = d0.mag2();
171 G4double sqrLen1 = d1.mag2(); 185 G4double sqrLen1 = d1.mag2();
172 location[0] = G4TwoVector(0.0,0.0); 186 location[0] = G4TwoVector(0.0,0.0);
173 location[1] = G4TwoVector(0.0,0.0); 187 location[1] = G4TwoVector(0.0,0.0);
174 188
175 if (sqrKross > DBL_EPSILON * DBL_EPSILON * s 189 if (sqrKross > DBL_EPSILON * DBL_EPSILON * sqrLen0 * sqrLen1)
176 { 190 {
177 // 191 //
178 // The line and line segment are not paral 192 // The line and line segment are not parallel. Determine if the intersection
179 // is in positive s where r=p0 + s*d0, and 193 // is in positive s where r=p0 + s*d0, and for 0<=t<=1 where r=p1 + t*d1.
180 // 194 //
181 G4double ss = cross(e,d1)/kross; 195 G4double ss = cross(e,d1)/kross;
182 if (ss < 0) return 0; // Intersect 196 if (ss < 0) return 0; // Intersection does not occur for positive ss
183 G4double t = cross(e,d0)/kross; 197 G4double t = cross(e,d0)/kross;
184 if (t < 0 || t > 1) return 0; // Intersect 198 if (t < 0 || t > 1) return 0; // Intersection does not occur on line-segment
185 // 199 //
186 // Intersection of lines is a single point 200 // Intersection of lines is a single point on the forward-propagating line
187 // defined by r=p0 + ss*d0, and the line s 201 // defined by r=p0 + ss*d0, and the line segment defined by r=p1 + t*d1.
188 // 202 //
189 location[0] = p0 + ss*d0; 203 location[0] = p0 + ss*d0;
190 return 1; 204 return 1;
191 } 205 }
192 // 206 //
193 // Line and line segment are parallel. Deter 207 // Line and line segment are parallel. Determine whether they overlap or not.
194 // 208 //
195 G4double sqrLenE = e.mag2(); 209 G4double sqrLenE = e.mag2();
196 kross = cross(e,d0); 210 kross = cross(e,d0);
197 sqrKross = kross * kross; 211 sqrKross = kross * kross;
198 if (sqrKross > DBL_EPSILON * DBL_EPSILON * s 212 if (sqrKross > DBL_EPSILON * DBL_EPSILON * sqrLen0 * sqrLenE)
199 { 213 {
200 return 0; //Lines are different. 214 return 0; //Lines are different.
201 } 215 }
202 // 216 //
203 // Lines are the same. Test for overlap. 217 // Lines are the same. Test for overlap.
204 // 218 //
205 G4double s0 = d0.dot(e)/sqrLen0; 219 G4double s0 = d0.dot(e)/sqrLen0;
206 G4double s1 = s0 + d0.dot(d1)/sqrLen0; 220 G4double s1 = s0 + d0.dot(d1)/sqrLen0;
207 G4double smin = 0.0; 221 G4double smin = 0.0;
208 G4double smax = 0.0; 222 G4double smax = 0.0;
209 223
210 if (s0 < s1) {smin = s0; smax = s1;} 224 if (s0 < s1) {smin = s0; smax = s1;}
211 else {smin = s1; smax = s0;} 225 else {smin = s1; smax = s0;}
212 226
213 if (smax < 0.0) return 0; 227 if (smax < 0.0) return 0;
214 else if (smin < 0.0) 228 else if (smin < 0.0)
215 { 229 {
216 location[0] = p0; 230 location[0] = p0;
217 location[1] = p0 + smax*d0; 231 location[1] = p0 + smax*d0;
218 return 2; 232 return 2;
219 } 233 }
220 else 234 else
221 { 235 {
222 location[0] = p0 + smin*d0; 236 location[0] = p0 + smin*d0;
223 location[1] = p0 + smax*d0; 237 location[1] = p0 + smax*d0;
224 return 2; 238 return 2;
225 } 239 }
226 } 240 }
227 241
228 ////////////////////////////////////////////// 242 ///////////////////////////////////////////////////////////////////////////////
229 // 243 //
230 // CrossProduct 244 // CrossProduct
231 // 245 //
232 // This is just a ficticious "cross-product" f 246 // This is just a ficticious "cross-product" function for two 2D vectors...
233 // "ficticious" because such an operation is n 247 // "ficticious" because such an operation is not relevant to 2D space compared
234 // with 3D space. 248 // with 3D space.
235 // 249 //
236 G4double G4TessellatedGeometryAlgorithms::cros << 250 G4double G4TessellatedGeometryAlgorithms::cross(const G4TwoVector &v1,
237 << 251 const G4TwoVector &v2)
238 { 252 {
239 return v1.x()*v2.y() - v1.y()*v2.x(); 253 return v1.x()*v2.y() - v1.y()*v2.x();
240 } 254 }
241 255