| 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. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // Implementation of methods for the class G4I 26 // Implementation of methods for the class G4IntersectionSolid
27 // 27 //
28 // 12.09.98 V.Grichine: first implementation 28 // 12.09.98 V.Grichine: first implementation
29 // ------------------------------------------- 29 // --------------------------------------------------------------------
30 30
31 #include <sstream> 31 #include <sstream>
32 32
33 #include "G4IntersectionSolid.hh" 33 #include "G4IntersectionSolid.hh"
34 34
35 #include "G4SystemOfUnits.hh" 35 #include "G4SystemOfUnits.hh"
36 #include "G4VoxelLimits.hh" 36 #include "G4VoxelLimits.hh"
37 #include "G4VPVParameterisation.hh" 37 #include "G4VPVParameterisation.hh"
38 38
39 #include "G4VGraphicsScene.hh" 39 #include "G4VGraphicsScene.hh"
40 #include "G4Polyhedron.hh" 40 #include "G4Polyhedron.hh"
41 #include "G4PolyhedronArbitrary.hh" 41 #include "G4PolyhedronArbitrary.hh"
42 #include "HepPolyhedronProcessor.h" 42 #include "HepPolyhedronProcessor.h"
43 43
44 ////////////////////////////////////////////// 44 //////////////////////////////////////////////////////////////////////////
45 // 45 //
46 // Transfer all data members to G4BooleanSolid 46 // Transfer all data members to G4BooleanSolid which is responsible
47 // for them. pName will be in turn sent to G4V 47 // for them. pName will be in turn sent to G4VSolid
48 // 48 //
49 49
50 G4IntersectionSolid::G4IntersectionSolid( cons 50 G4IntersectionSolid::G4IntersectionSolid( const G4String& pName,
51 51 G4VSolid* pSolidA ,
52 52 G4VSolid* pSolidB )
53 : G4BooleanSolid(pName,pSolidA,pSolidB) 53 : G4BooleanSolid(pName,pSolidA,pSolidB)
54 { 54 {
55 } 55 }
56 56
57 ////////////////////////////////////////////// 57 //////////////////////////////////////////////////////////////////////////
58 // 58 //
59 59
60 G4IntersectionSolid::G4IntersectionSolid( cons 60 G4IntersectionSolid::G4IntersectionSolid( const G4String& pName,
61 61 G4VSolid* pSolidA,
62 62 G4VSolid* pSolidB,
63 63 G4RotationMatrix* rotMatrix,
64 cons 64 const G4ThreeVector& transVector )
65 : G4BooleanSolid(pName,pSolidA,pSolidB,rotMa 65 : G4BooleanSolid(pName,pSolidA,pSolidB,rotMatrix,transVector)
66 { 66 {
67 } 67 }
68 68
69 ////////////////////////////////////////////// 69 //////////////////////////////////////////////////////////////////////////
70 // 70 //
71 // 71 //
72 72
73 G4IntersectionSolid::G4IntersectionSolid( cons 73 G4IntersectionSolid::G4IntersectionSolid( const G4String& pName,
74 74 G4VSolid* pSolidA,
75 75 G4VSolid* pSolidB,
76 cons 76 const G4Transform3D& transform )
77 : G4BooleanSolid(pName,pSolidA,pSolidB,trans 77 : G4BooleanSolid(pName,pSolidA,pSolidB,transform)
78 { 78 {
79 } 79 }
80 80
81 ////////////////////////////////////////////// 81 //////////////////////////////////////////////////////////////////////////
82 // 82 //
83 // Fake default constructor - sets only member 83 // Fake default constructor - sets only member data and allocates memory
84 // for usage restri 84 // for usage restricted to object persistency.
85 85
86 G4IntersectionSolid::G4IntersectionSolid( __vo 86 G4IntersectionSolid::G4IntersectionSolid( __void__& a )
87 : G4BooleanSolid(a) 87 : G4BooleanSolid(a)
88 { 88 {
89 } 89 }
90 90
91 ////////////////////////////////////////////// 91 //////////////////////////////////////////////////////////////////////////
92 // 92 //
93 // 93 //
94 94
95 G4IntersectionSolid::~G4IntersectionSolid() = 95 G4IntersectionSolid::~G4IntersectionSolid() = default;
96 96
97 ////////////////////////////////////////////// 97 //////////////////////////////////////////////////////////////////////////
98 // 98 //
99 // Copy constructor 99 // Copy constructor
100 100
101 G4IntersectionSolid::G4IntersectionSolid(const 101 G4IntersectionSolid::G4IntersectionSolid(const G4IntersectionSolid&) = default;
102 102
103 ////////////////////////////////////////////// 103 //////////////////////////////////////////////////////////////////////////
104 // 104 //
105 // Assignment operator 105 // Assignment operator
106 106
107 G4IntersectionSolid& 107 G4IntersectionSolid&
108 G4IntersectionSolid::operator = (const G4Inter 108 G4IntersectionSolid::operator = (const G4IntersectionSolid& rhs)
109 { 109 {
110 // Check assignment to self 110 // Check assignment to self
111 // 111 //
112 if (this == &rhs) { return *this; } 112 if (this == &rhs) { return *this; }
113 113
114 // Copy base class data 114 // Copy base class data
115 // 115 //
116 G4BooleanSolid::operator=(rhs); 116 G4BooleanSolid::operator=(rhs);
117 117
118 return *this; 118 return *this;
119 } 119 }
120 120
121 ////////////////////////////////////////////// 121 //////////////////////////////////////////////////////////////////////////
122 // 122 //
123 // Get bounding box 123 // Get bounding box
124 124
125 void 125 void
126 G4IntersectionSolid::BoundingLimits(G4ThreeVec 126 G4IntersectionSolid::BoundingLimits(G4ThreeVector& pMin,
127 G4ThreeVec 127 G4ThreeVector& pMax) const
128 { 128 {
129 G4ThreeVector minA,maxA, minB,maxB; 129 G4ThreeVector minA,maxA, minB,maxB;
130 fPtrSolidA->BoundingLimits(minA,maxA); 130 fPtrSolidA->BoundingLimits(minA,maxA);
131 fPtrSolidB->BoundingLimits(minB,maxB); 131 fPtrSolidB->BoundingLimits(minB,maxB);
132 132
133 pMin.set(std::max(minA.x(),minB.x()), 133 pMin.set(std::max(minA.x(),minB.x()),
134 std::max(minA.y(),minB.y()), 134 std::max(minA.y(),minB.y()),
135 std::max(minA.z(),minB.z())); 135 std::max(minA.z(),minB.z()));
136 136
137 pMax.set(std::min(maxA.x(),maxB.x()), 137 pMax.set(std::min(maxA.x(),maxB.x()),
138 std::min(maxA.y(),maxB.y()), 138 std::min(maxA.y(),maxB.y()),
139 std::min(maxA.z(),maxB.z())); 139 std::min(maxA.z(),maxB.z()));
140 140
141 // Check correctness of the bounding box 141 // Check correctness of the bounding box
142 // 142 //
143 if (pMin.x() >= pMax.x() || pMin.y() >= pMax 143 if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
144 { 144 {
145 std::ostringstream message; 145 std::ostringstream message;
146 message << "Bad bounding box (min >= max) 146 message << "Bad bounding box (min >= max) for solid: "
147 << GetName() << " !" 147 << GetName() << " !"
148 << "\npMin = " << pMin 148 << "\npMin = " << pMin
149 << "\npMax = " << pMax; 149 << "\npMax = " << pMax;
150 G4Exception("G4IntersectionSolid::Bounding 150 G4Exception("G4IntersectionSolid::BoundingLimits()", "GeomMgt0001",
151 JustWarning, message); 151 JustWarning, message);
152 DumpInfo(); 152 DumpInfo();
153 } 153 }
154 } 154 }
155 155
156 ////////////////////////////////////////////// 156 //////////////////////////////////////////////////////////////////////////
157 // 157 //
158 // Calculate extent under transform and specif 158 // Calculate extent under transform and specified limit
159 159
160 G4bool 160 G4bool
161 G4IntersectionSolid::CalculateExtent(const EAx 161 G4IntersectionSolid::CalculateExtent(const EAxis pAxis,
162 const G4V 162 const G4VoxelLimits& pVoxelLimit,
163 const G4A 163 const G4AffineTransform& pTransform,
164 G4d 164 G4double& pMin,
165 G4d 165 G4double& pMax) const
166 { 166 {
167 G4bool retA, retB, out; 167 G4bool retA, retB, out;
168 G4double minA, minB, maxA, maxB; 168 G4double minA, minB, maxA, maxB;
169 169
170 retA = fPtrSolidA 170 retA = fPtrSolidA
171 ->CalculateExtent( pAxis, pVoxelLimi 171 ->CalculateExtent( pAxis, pVoxelLimit, pTransform, minA, maxA);
172 retB = fPtrSolidB 172 retB = fPtrSolidB
173 ->CalculateExtent( pAxis, pVoxelLimi 173 ->CalculateExtent( pAxis, pVoxelLimit, pTransform, minB, maxB);
174 174
175 if( retA && retB ) 175 if( retA && retB )
176 { 176 {
177 pMin = std::max( minA, minB ); 177 pMin = std::max( minA, minB );
178 pMax = std::min( maxA, maxB ); 178 pMax = std::min( maxA, maxB );
179 out = (pMax > pMin); // true; 179 out = (pMax > pMin); // true;
180 } 180 }
181 else 181 else
182 { 182 {
183 out = false; 183 out = false;
184 } 184 }
185 185
186 return out; // It exists in this slice only 186 return out; // It exists in this slice only if both exist in it.
187 } 187 }
188 188
189 ////////////////////////////////////////////// 189 //////////////////////////////////////////////////////////////////////////
190 // 190 //
191 // Touching ? Empty intersection ? 191 // Touching ? Empty intersection ?
192 192
193 EInside G4IntersectionSolid::Inside(const G4Th 193 EInside G4IntersectionSolid::Inside(const G4ThreeVector& p) const
194 { 194 {
195 EInside positionA = fPtrSolidA->Inside(p); 195 EInside positionA = fPtrSolidA->Inside(p);
196 if(positionA == kOutside) return positionA; 196 if(positionA == kOutside) return positionA; // outside A
197 197
198 EInside positionB = fPtrSolidB->Inside(p); 198 EInside positionB = fPtrSolidB->Inside(p);
199 if(positionA == kInside) return positionB; 199 if(positionA == kInside) return positionB;
200 200
201 if(positionB == kOutside) return positionB; 201 if(positionB == kOutside) return positionB; // outside B
202 return kSurface; 202 return kSurface; // surface A & B
203 } 203 }
204 204
205 ////////////////////////////////////////////// 205 //////////////////////////////////////////////////////////////////////////
206 // 206 //
207 207
208 G4ThreeVector 208 G4ThreeVector
209 G4IntersectionSolid::SurfaceNormal( const G4Th 209 G4IntersectionSolid::SurfaceNormal( const G4ThreeVector& p ) const
210 { 210 {
211 G4ThreeVector normal; 211 G4ThreeVector normal;
212 EInside insideA, insideB; 212 EInside insideA, insideB;
213 213
214 insideA = fPtrSolidA->Inside(p); 214 insideA = fPtrSolidA->Inside(p);
215 insideB = fPtrSolidB->Inside(p); 215 insideB = fPtrSolidB->Inside(p);
216 216
217 #ifdef G4BOOLDEBUG 217 #ifdef G4BOOLDEBUG
218 if( (insideA == kOutside) || (insideB == kOu 218 if( (insideA == kOutside) || (insideB == kOutside) )
219 { 219 {
220 G4cout << "WARNING - Invalid call in " 220 G4cout << "WARNING - Invalid call in "
221 << "G4IntersectionSolid::SurfaceNor 221 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
222 << " Point p is outside !" << G4en 222 << " Point p is outside !" << G4endl;
223 G4cout << " p = " << p << G4endl; 223 G4cout << " p = " << p << G4endl;
224 G4cerr << "WARNING - Invalid call in " 224 G4cerr << "WARNING - Invalid call in "
225 << "G4IntersectionSolid::SurfaceNor 225 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
226 << " Point p is outside !" << G4en 226 << " Point p is outside !" << G4endl;
227 G4cerr << " p = " << p << G4endl; 227 G4cerr << " p = " << p << G4endl;
228 } 228 }
229 #endif 229 #endif
230 230
231 // On the surface of both is difficult ... t 231 // On the surface of both is difficult ... treat it like on A now!
232 // 232 //
233 if( insideA == kSurface ) 233 if( insideA == kSurface )
234 { 234 {
235 normal = fPtrSolidA->SurfaceNormal(p) ; 235 normal = fPtrSolidA->SurfaceNormal(p) ;
236 } 236 }
237 else if( insideB == kSurface ) 237 else if( insideB == kSurface )
238 { 238 {
239 normal = fPtrSolidB->SurfaceNormal(p) ; 239 normal = fPtrSolidB->SurfaceNormal(p) ;
240 } 240 }
241 else // We are on neither surface, so we sh 241 else // We are on neither surface, so we should generate an exception
242 { 242 {
243 if(fPtrSolidA->DistanceToOut(p) <= fPtrSol 243 if(fPtrSolidA->DistanceToOut(p) <= fPtrSolidB->DistanceToOut(p) )
244 { 244 {
245 normal= fPtrSolidA->SurfaceNormal(p) ; 245 normal= fPtrSolidA->SurfaceNormal(p) ;
246 } 246 }
247 else 247 else
248 { 248 {
249 normal= fPtrSolidB->SurfaceNormal(p) ; 249 normal= fPtrSolidB->SurfaceNormal(p) ;
250 } 250 }
251 #ifdef G4BOOLDEBUG 251 #ifdef G4BOOLDEBUG
252 G4cout << "WARNING - Invalid call in " 252 G4cout << "WARNING - Invalid call in "
253 << "G4IntersectionSolid::SurfaceNor 253 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
254 << " Point p is out of surface !" 254 << " Point p is out of surface !" << G4endl;
255 G4cout << " p = " << p << G4endl; 255 G4cout << " p = " << p << G4endl;
256 G4cerr << "WARNING - Invalid call in " 256 G4cerr << "WARNING - Invalid call in "
257 << "G4IntersectionSolid::SurfaceNor 257 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
258 << " Point p is out of surface !" 258 << " Point p is out of surface !" << G4endl;
259 G4cerr << " p = " << p << G4endl; 259 G4cerr << " p = " << p << G4endl;
260 #endif 260 #endif
261 } 261 }
262 262
263 return normal; 263 return normal;
264 } 264 }
265 265
266 ////////////////////////////////////////////// 266 //////////////////////////////////////////////////////////////////////////
267 // 267 //
268 // The same algorithm as in DistanceToIn(p) 268 // The same algorithm as in DistanceToIn(p)
269 269
270 G4double 270 G4double
271 G4IntersectionSolid::DistanceToIn( const G4Thr 271 G4IntersectionSolid::DistanceToIn( const G4ThreeVector& p,
272 const G4Thr 272 const G4ThreeVector& v ) const
273 { 273 {
274 G4double dist = 0.0; 274 G4double dist = 0.0;
275 if( Inside(p) == kInside ) 275 if( Inside(p) == kInside )
276 { 276 {
277 #ifdef G4BOOLDEBUG 277 #ifdef G4BOOLDEBUG
278 G4cout << "WARNING - Invalid call in " 278 G4cout << "WARNING - Invalid call in "
279 << "G4IntersectionSolid::DistanceTo 279 << "G4IntersectionSolid::DistanceToIn(p,v)" << G4endl
280 << " Point p is inside !" << G4end 280 << " Point p is inside !" << G4endl;
281 G4cout << " p = " << p << G4endl; 281 G4cout << " p = " << p << G4endl;
282 G4cout << " v = " << v << G4endl; 282 G4cout << " v = " << v << G4endl;
283 G4cerr << "WARNING - Invalid call in " 283 G4cerr << "WARNING - Invalid call in "
284 << "G4IntersectionSolid::DistanceTo 284 << "G4IntersectionSolid::DistanceToIn(p,v)" << G4endl
285 << " Point p is inside !" << G4end 285 << " Point p is inside !" << G4endl;
286 G4cerr << " p = " << p << G4endl; 286 G4cerr << " p = " << p << G4endl;
287 G4cerr << " v = " << v << G4endl; 287 G4cerr << " v = " << v << G4endl;
288 #endif 288 #endif
289 } 289 }
290 else // if( Inside(p) == kSurface ) 290 else // if( Inside(p) == kSurface )
291 { 291 {
292 EInside wA = fPtrSolidA->Inside(p); 292 EInside wA = fPtrSolidA->Inside(p);
293 EInside wB = fPtrSolidB->Inside(p); 293 EInside wB = fPtrSolidB->Inside(p);
294 294
295 G4ThreeVector pA = p, pB = p; 295 G4ThreeVector pA = p, pB = p;
296 G4double dA = 0., dA1=0., dA2=0.; 296 G4double dA = 0., dA1=0., dA2=0.;
297 G4double dB = 0., dB1=0., dB2=0.; 297 G4double dB = 0., dB1=0., dB2=0.;
298 G4bool doA = true, doB = true; 298 G4bool doA = true, doB = true;
299 299
300 static const std::size_t max_trials=10000; 300 static const std::size_t max_trials=10000;
301 for (std::size_t trial=0; trial<max_trials 301 for (std::size_t trial=0; trial<max_trials; ++trial)
302 { 302 {
303 if(doA) 303 if(doA)
304 { 304 {
305 // find next valid range for A 305 // find next valid range for A
306 306
307 dA1 = 0.; 307 dA1 = 0.;
308 308
309 if( wA != kInside ) 309 if( wA != kInside )
310 { 310 {
311 dA1 = fPtrSolidA->DistanceToIn(pA, v 311 dA1 = fPtrSolidA->DistanceToIn(pA, v);
312 312
313 if( dA1 == kInfinity ) return kInf 313 if( dA1 == kInfinity ) return kInfinity;
314 314
315 pA += dA1*v; 315 pA += dA1*v;
316 } 316 }
317 dA2 = dA1 + fPtrSolidA->DistanceToOut( 317 dA2 = dA1 + fPtrSolidA->DistanceToOut(pA, v);
318 } 318 }
319 dA1 += dA; 319 dA1 += dA;
320 dA2 += dA; 320 dA2 += dA;
321 321
322 if(doB) 322 if(doB)
323 { 323 {
324 // find next valid range for B 324 // find next valid range for B
325 325
326 dB1 = 0.; 326 dB1 = 0.;
327 if(wB != kInside) 327 if(wB != kInside)
328 { 328 {
329 dB1 = fPtrSolidB->DistanceToIn(pB, v 329 dB1 = fPtrSolidB->DistanceToIn(pB, v);
330 330
331 if(dB1 == kInfinity) return kInfin 331 if(dB1 == kInfinity) return kInfinity;
332 332
333 pB += dB1*v; 333 pB += dB1*v;
334 } 334 }
335 dB2 = dB1 + fPtrSolidB->DistanceToOut( 335 dB2 = dB1 + fPtrSolidB->DistanceToOut(pB, v);
336 } 336 }
337 dB1 += dB; 337 dB1 += dB;
338 dB2 += dB; 338 dB2 += dB;
339 339
340 // check if they overlap 340 // check if they overlap
341 341
342 if( dA1 < dB1 ) 342 if( dA1 < dB1 )
343 { 343 {
344 if( dB1 < dA2 ) return dB1; 344 if( dB1 < dA2 ) return dB1;
345 345
346 dA = dA2; 346 dA = dA2;
347 pA = p + dA*v; // continue from her 347 pA = p + dA*v; // continue from here
348 wA = kSurface; 348 wA = kSurface;
349 doA = true; 349 doA = true;
350 doB = false; 350 doB = false;
351 } 351 }
352 else 352 else
353 { 353 {
354 if( dA1 < dB2 ) return dA1; 354 if( dA1 < dB2 ) return dA1;
355 355
356 dB = dB2; 356 dB = dB2;
357 pB = p + dB*v; // continue from her 357 pB = p + dB*v; // continue from here
358 wB = kSurface; 358 wB = kSurface;
359 doB = true; 359 doB = true;
360 doA = false; 360 doA = false;
361 } 361 }
362 } 362 }
363 } 363 }
364 #ifdef G4BOOLDEBUG 364 #ifdef G4BOOLDEBUG
365 G4Exception("G4IntersectionSolid::DistanceTo 365 G4Exception("G4IntersectionSolid::DistanceToIn(p,v)",
366 "GeomSolids0001", JustWarning, 366 "GeomSolids0001", JustWarning,
367 "Reached maximum number of itera 367 "Reached maximum number of iterations! Returning zero.");
368 #endif 368 #endif
369 return dist ; 369 return dist ;
370 } 370 }
371 371
372 ////////////////////////////////////////////// 372 //////////////////////////////////////////////////////////////////////////
373 // 373 //
374 // Approximate nearest distance from the point 374 // Approximate nearest distance from the point p to the intersection of
375 // two solids 375 // two solids
376 376
377 G4double 377 G4double
378 G4IntersectionSolid::DistanceToIn( const G4Thr 378 G4IntersectionSolid::DistanceToIn( const G4ThreeVector& p) const
379 { 379 {
380 #ifdef G4BOOLDEBUG 380 #ifdef G4BOOLDEBUG
381 if( Inside(p) == kInside ) 381 if( Inside(p) == kInside )
382 { 382 {
383 G4cout << "WARNING - Invalid call in " 383 G4cout << "WARNING - Invalid call in "
384 << "G4IntersectionSolid::DistanceTo 384 << "G4IntersectionSolid::DistanceToIn(p)" << G4endl
385 << " Point p is inside !" << G4end 385 << " Point p is inside !" << G4endl;
386 G4cout << " p = " << p << G4endl; 386 G4cout << " p = " << p << G4endl;
387 G4cerr << "WARNING - Invalid call in " 387 G4cerr << "WARNING - Invalid call in "
388 << "G4IntersectionSolid::DistanceTo 388 << "G4IntersectionSolid::DistanceToIn(p)" << G4endl
389 << " Point p is inside !" << G4end 389 << " Point p is inside !" << G4endl;
390 G4cerr << " p = " << p << G4endl; 390 G4cerr << " p = " << p << G4endl;
391 } 391 }
392 #endif 392 #endif
393 EInside sideA = fPtrSolidA->Inside(p) ; 393 EInside sideA = fPtrSolidA->Inside(p) ;
394 EInside sideB = fPtrSolidB->Inside(p) ; 394 EInside sideB = fPtrSolidB->Inside(p) ;
395 G4double dist=0.0 ; 395 G4double dist=0.0 ;
396 396
397 if( sideA != kInside && sideB != kOutside ) 397 if( sideA != kInside && sideB != kOutside )
398 { 398 {
399 dist = fPtrSolidA->DistanceToIn(p) ; 399 dist = fPtrSolidA->DistanceToIn(p) ;
400 } 400 }
401 else 401 else
402 { 402 {
403 if( sideB != kInside && sideA != kOutside 403 if( sideB != kInside && sideA != kOutside )
404 { 404 {
405 dist = fPtrSolidB->DistanceToIn(p) ; 405 dist = fPtrSolidB->DistanceToIn(p) ;
406 } 406 }
407 else 407 else
408 { 408 {
409 dist = std::min(fPtrSolidA->DistanceToI 409 dist = std::min(fPtrSolidA->DistanceToIn(p),
410 fPtrSolidB->DistanceToI 410 fPtrSolidB->DistanceToIn(p) ) ;
411 } 411 }
412 } 412 }
413 return dist ; 413 return dist ;
414 } 414 }
415 415
416 ////////////////////////////////////////////// 416 //////////////////////////////////////////////////////////////////////////
417 // 417 //
418 // The same algorithm as DistanceToOut(p) 418 // The same algorithm as DistanceToOut(p)
419 419
420 G4double 420 G4double
421 G4IntersectionSolid::DistanceToOut( const G4Th 421 G4IntersectionSolid::DistanceToOut( const G4ThreeVector& p,
422 const G4Th 422 const G4ThreeVector& v,
423 const G4bo 423 const G4bool calcNorm,
424 G4bo 424 G4bool* validNorm,
425 G4Th 425 G4ThreeVector* n ) const
426 { 426 {
427 G4bool validNormA, validNormB; 427 G4bool validNormA, validNormB;
428 G4ThreeVector nA, nB; 428 G4ThreeVector nA, nB;
429 429
430 #ifdef G4BOOLDEBUG 430 #ifdef G4BOOLDEBUG
431 if( Inside(p) == kOutside ) 431 if( Inside(p) == kOutside )
432 { 432 {
433 G4cout << "Position:" << G4endl << G4endl 433 G4cout << "Position:" << G4endl << G4endl;
434 G4cout << "p.x() = " << p.x()/mm << " mm 434 G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl;
435 G4cout << "p.y() = " << p.y()/mm << " mm 435 G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl;
436 G4cout << "p.z() = " << p.z()/mm << " mm 436 G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl;
437 G4cout << "Direction:" << G4endl << G4endl 437 G4cout << "Direction:" << G4endl << G4endl;
438 G4cout << "v.x() = " << v.x() << G4endl; 438 G4cout << "v.x() = " << v.x() << G4endl;
439 G4cout << "v.y() = " << v.y() << G4endl; 439 G4cout << "v.y() = " << v.y() << G4endl;
440 G4cout << "v.z() = " << v.z() << G4endl 440 G4cout << "v.z() = " << v.z() << G4endl << G4endl;
441 G4cout << "WARNING - Invalid call in " 441 G4cout << "WARNING - Invalid call in "
442 << "G4IntersectionSolid::DistanceTo 442 << "G4IntersectionSolid::DistanceToOut(p,v)" << G4endl
443 << " Point p is outside !" << G4en 443 << " Point p is outside !" << G4endl;
444 G4cout << " p = " << p << G4endl; 444 G4cout << " p = " << p << G4endl;
445 G4cout << " v = " << v << G4endl; 445 G4cout << " v = " << v << G4endl;
446 G4cerr << "WARNING - Invalid call in " 446 G4cerr << "WARNING - Invalid call in "
447 << "G4IntersectionSolid::DistanceTo 447 << "G4IntersectionSolid::DistanceToOut(p,v)" << G4endl
448 << " Point p is outside !" << G4en 448 << " Point p is outside !" << G4endl;
449 G4cerr << " p = " << p << G4endl; 449 G4cerr << " p = " << p << G4endl;
450 G4cerr << " v = " << v << G4endl; 450 G4cerr << " v = " << v << G4endl;
451 } 451 }
452 #endif 452 #endif
453 G4double distA = fPtrSolidA->DistanceToOut(p 453 G4double distA = fPtrSolidA->DistanceToOut(p,v,calcNorm,&validNormA,&nA) ;
454 G4double distB = fPtrSolidB->DistanceToOut(p 454 G4double distB = fPtrSolidB->DistanceToOut(p,v,calcNorm,&validNormB,&nB) ;
455 455
456 G4double dist = std::min(distA,distB) ; 456 G4double dist = std::min(distA,distB) ;
457 457
458 if( calcNorm ) 458 if( calcNorm )
459 { 459 {
460 if ( distA < distB ) 460 if ( distA < distB )
461 { 461 {
462 *validNorm = validNormA; 462 *validNorm = validNormA;
463 *n = nA; 463 *n = nA;
464 } 464 }
465 else 465 else
466 { 466 {
467 *validNorm = validNormB; 467 *validNorm = validNormB;
468 *n = nB; 468 *n = nB;
469 } 469 }
470 } 470 }
471 471
472 return dist ; 472 return dist ;
473 } 473 }
474 474
475 ////////////////////////////////////////////// 475 //////////////////////////////////////////////////////////////////////////
476 // 476 //
477 // Inverted algorithm of DistanceToIn(p) 477 // Inverted algorithm of DistanceToIn(p)
478 478
479 G4double 479 G4double
480 G4IntersectionSolid::DistanceToOut( const G4Th 480 G4IntersectionSolid::DistanceToOut( const G4ThreeVector& p ) const
481 { 481 {
482 #ifdef G4BOOLDEBUG 482 #ifdef G4BOOLDEBUG
483 if( Inside(p) == kOutside ) 483 if( Inside(p) == kOutside )
484 { 484 {
485 G4cout << "WARNING - Invalid call in " 485 G4cout << "WARNING - Invalid call in "
486 << "G4IntersectionSolid::DistanceTo 486 << "G4IntersectionSolid::DistanceToOut(p)" << G4endl
487 << " Point p is outside !" << G4en 487 << " Point p is outside !" << G4endl;
488 G4cout << " p = " << p << G4endl; 488 G4cout << " p = " << p << G4endl;
489 G4cerr << "WARNING - Invalid call in " 489 G4cerr << "WARNING - Invalid call in "
490 << "G4IntersectionSolid::DistanceTo 490 << "G4IntersectionSolid::DistanceToOut(p)" << G4endl
491 << " Point p is outside !" << G4en 491 << " Point p is outside !" << G4endl;
492 G4cerr << " p = " << p << G4endl; 492 G4cerr << " p = " << p << G4endl;
493 } 493 }
494 #endif 494 #endif
495 495
496 return std::min(fPtrSolidA->DistanceToOut(p) 496 return std::min(fPtrSolidA->DistanceToOut(p),
497 fPtrSolidB->DistanceToOut(p) 497 fPtrSolidB->DistanceToOut(p) ) ;
498 498
499 } 499 }
500 500
501 ////////////////////////////////////////////// 501 //////////////////////////////////////////////////////////////////////////
502 // 502 //
503 // ComputeDimensions 503 // ComputeDimensions
504 504
505 void 505 void
506 G4IntersectionSolid::ComputeDimensions( G4VPVP 506 G4IntersectionSolid::ComputeDimensions( G4VPVParameterisation*,
507 const 507 const G4int,
508 const 508 const G4VPhysicalVolume* )
509 { 509 {
510 } 510 }
511 511
512 ////////////////////////////////////////////// 512 //////////////////////////////////////////////////////////////////////////
513 // 513 //
514 // GetEntityType 514 // GetEntityType
515 515
516 G4GeometryType G4IntersectionSolid::GetEntityT 516 G4GeometryType G4IntersectionSolid::GetEntityType() const
517 { 517 {
518 return {"G4IntersectionSolid"}; 518 return {"G4IntersectionSolid"};
519 } 519 }
520 520
521 ////////////////////////////////////////////// 521 //////////////////////////////////////////////////////////////////////////
522 // 522 //
523 // Make a clone of the object 523 // Make a clone of the object
524 524
525 G4VSolid* G4IntersectionSolid::Clone() const 525 G4VSolid* G4IntersectionSolid::Clone() const
526 { 526 {
527 return new G4IntersectionSolid(*this); 527 return new G4IntersectionSolid(*this);
528 } 528 }
529 529
530 ////////////////////////////////////////////// 530 //////////////////////////////////////////////////////////////////////////
531 // 531 //
532 // DescribeYourselfTo 532 // DescribeYourselfTo
533 533
534 void 534 void
535 G4IntersectionSolid::DescribeYourselfTo ( G4VG 535 G4IntersectionSolid::DescribeYourselfTo ( G4VGraphicsScene& scene ) const
536 { 536 {
537 scene.AddSolid (*this); 537 scene.AddSolid (*this);
538 } 538 }
539 539
540 ////////////////////////////////////////////// 540 //////////////////////////////////////////////////////////////////////////
541 // 541 //
542 // CreatePolyhedron 542 // CreatePolyhedron
543 543
544 G4Polyhedron* 544 G4Polyhedron*
545 G4IntersectionSolid::CreatePolyhedron () const 545 G4IntersectionSolid::CreatePolyhedron () const
546 { 546 {
547 if (fExternalBoolProcessor == nullptr) 547 if (fExternalBoolProcessor == nullptr)
548 { 548 {
549 HepPolyhedronProcessor processor; 549 HepPolyhedronProcessor processor;
550 // Stack components and components of comp 550 // Stack components and components of components recursively
551 // See G4BooleanSolid::StackPolyhedron 551 // See G4BooleanSolid::StackPolyhedron
552 G4Polyhedron* top = StackPolyhedron(proces 552 G4Polyhedron* top = StackPolyhedron(processor, this);
553 auto result = new G4Polyhedron(*top); 553 auto result = new G4Polyhedron(*top);
554 if (processor.execute(*result)) 554 if (processor.execute(*result))
555 { 555 {
556 return result; 556 return result;
557 } 557 }
558 else 558 else
559 { 559 {
560 return nullptr; 560 return nullptr;
561 } 561 }
562 } 562 }
563 else 563 else
564 { 564 {
565 return fExternalBoolProcessor->Process(thi << 565 return fExternalBoolProcessor
>> 566 ->Intersection(GetConstituentSolid(0)->GetPolyhedron(),
>> 567 GetConstituentSolid(1)->GetPolyhedron());
566 } 568 }
567 } 569 }
568 570