| 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 //
>> 27 // $Id: G4NavigationLogger.cc 93289 2015-10-15 10:01:15Z gcosmo $
>> 28 //
>> 29 //
26 // class G4NavigationLogger Implementation 30 // class G4NavigationLogger Implementation
27 // 31 //
28 // Author: G.Cosmo, 2010 32 // Author: G.Cosmo, 2010
>> 33 //
29 // ------------------------------------------- 34 // --------------------------------------------------------------------
30 35
31 #include <iomanip> 36 #include <iomanip>
32 #include <CLHEP/Units/SystemOfUnits.h> 37 #include <CLHEP/Units/SystemOfUnits.h>
33 38
34 #include "G4NavigationLogger.hh" 39 #include "G4NavigationLogger.hh"
35 #include "G4GeometryTolerance.hh" 40 #include "G4GeometryTolerance.hh"
36 41
>> 42 // const double millimeter = CLHEP::Units::millimeter;
37 using CLHEP::millimeter; 43 using CLHEP::millimeter;
38 44
39 G4NavigationLogger::G4NavigationLogger(const G 45 G4NavigationLogger::G4NavigationLogger(const G4String& id)
40 : fId(id) << 46 : fId(id), fVerbose(0), fMinTriggerDistance( DBL_MAX ), fReportSoftWarnings( false )
41 { 47 {
42 } 48 }
43 49
44 G4NavigationLogger::~G4NavigationLogger() = de << 50 G4NavigationLogger::~G4NavigationLogger()
>> 51 {
>> 52 }
45 53
46 // ******************************************* 54 // ********************************************************************
47 // PreComputeStepLog 55 // PreComputeStepLog
48 // ******************************************* 56 // ********************************************************************
49 // 57 //
50 void 58 void
51 G4NavigationLogger::PreComputeStepLog(const G4 59 G4NavigationLogger::PreComputeStepLog(const G4VPhysicalVolume* motherPhysical,
52 G4 60 G4double motherSafety,
53 const G4 61 const G4ThreeVector& localPoint) const
54 { 62 {
55 G4VSolid* motherSolid = motherPhysical->GetL 63 G4VSolid* motherSolid = motherPhysical->GetLogicalVolume()->GetSolid();
56 G4String fType = fId + "::ComputeStep()"; 64 G4String fType = fId + "::ComputeStep()";
>> 65 // const double millimeter = CLHEP::millimeter;
57 66
58 if ( fVerbose == 1 || fVerbose > 4 ) 67 if ( fVerbose == 1 || fVerbose > 4 )
59 { 68 {
60 G4cout << "*************** " << fType << " 69 G4cout << "*************** " << fType << " *****************" << G4endl
61 << " VolType " 70 << " VolType "
62 << std::setw(15) << "Safety/mm" << 71 << std::setw(15) << "Safety/mm" << " "
63 << std::setw(15) << "Distance/mm" < 72 << std::setw(15) << "Distance/mm" << " "
64 << std::setw(52) << "Position (loca 73 << std::setw(52) << "Position (local coordinates)"
65 << " - Solid" << G4endl; 74 << " - Solid" << G4endl;
66 G4cout << " Mother " 75 G4cout << " Mother "
67 << std::setw(15) << motherSafety / 76 << std::setw(15) << motherSafety / millimeter << " "
68 << std::setw(15) << "N/C" << 77 << std::setw(15) << "N/C" << " " << localPoint << " - "
69 << motherSolid->GetEntityType() << 78 << motherSolid->GetEntityType() << ": " << motherSolid->GetName()
70 << G4endl; 79 << G4endl;
71 } 80 }
72 if ( motherSafety < 0.0 ) 81 if ( motherSafety < 0.0 )
73 { 82 {
74 std::ostringstream message; 83 std::ostringstream message;
75 message << "Negative Safety In Voxel Navig 84 message << "Negative Safety In Voxel Navigation !" << G4endl
76 << " Current solid " << mot 85 << " Current solid " << motherSolid->GetName()
77 << " gave negative safety: " << mo 86 << " gave negative safety: " << motherSafety / millimeter << G4endl
78 << " for the current (local 87 << " for the current (local) point " << localPoint;
79 message << " Solid info: " << *motherSolid 88 message << " Solid info: " << *motherSolid << G4endl;
80 G4Exception(fType, "GeomNav0003", FatalExc 89 G4Exception(fType, "GeomNav0003", FatalException, message);
81 } 90 }
82 if( motherSolid->Inside(localPoint)==kOutsid 91 if( motherSolid->Inside(localPoint)==kOutside )
83 { 92 {
84 std::ostringstream message; 93 std::ostringstream message;
85 message << "Point is outside Current Volum 94 message << "Point is outside Current Volume - " << G4endl
86 << " Point " << localPoin 95 << " Point " << localPoint / millimeter
87 << " is outside current volume '" 96 << " is outside current volume '" << motherPhysical->GetName()
88 << "'" << G4endl; 97 << "'" << G4endl;
89 G4double estDistToSolid= motherSolid->Dist 98 G4double estDistToSolid= motherSolid->DistanceToIn(localPoint);
90 message << " Estimated isotropic 99 message << " Estimated isotropic distance to solid (distToIn)= "
91 << estDistToSolid << G4endl; 100 << estDistToSolid << G4endl;
92 if( estDistToSolid > 100.0 * motherSolid-> 101 if( estDistToSolid > 100.0 * motherSolid->GetTolerance() )
93 { 102 {
94 message << " Solid info: " << *motherSol 103 message << " Solid info: " << *motherSolid << G4endl;
95 G4Exception(fType, "GeomNav0003", JustWa 104 G4Exception(fType, "GeomNav0003", JustWarning, message,
96 "Point is far outside Curren 105 "Point is far outside Current Volume !" );
97 } 106 }
98 else 107 else
99 { 108 {
100 G4Exception(fType, "GeomNav1001", JustWa 109 G4Exception(fType, "GeomNav1001", JustWarning, message,
101 "Point is a little outside C 110 "Point is a little outside Current Volume.");
102 } 111 }
103 } 112 }
104 113
105 // Verification / verbosity 114 // Verification / verbosity
106 // 115 //
107 if ( fVerbose > 1 ) 116 if ( fVerbose > 1 )
108 { 117 {
109 static const G4int precVerf = 16; // Prec << 118 static const G4int precVerf= 16; // Precision
110 G4long oldprec = G4cout.precision(precVerf << 119 G4int oldprec = G4cout.precision(precVerf);
111 G4cout << " - Information on mother / key 120 G4cout << " - Information on mother / key daughters ..." << G4endl;
112 G4cout << " Type " << std::setw(12) << 121 G4cout << " Type " << std::setw(12) << "Solid-Name" << " "
113 << std::setw(3*(6+precVerf)) << 122 << std::setw(3*(6+precVerf)) << " local point" << " "
114 << std::setw(4+precVerf) << 123 << std::setw(4+precVerf) << "solid-Safety" << " "
115 << std::setw(4+precVerf) << 124 << std::setw(4+precVerf) << "solid-Step" << " "
116 << std::setw(17) << 125 << std::setw(17) << "distance Method "
117 << std::setw(3*(6+precVerf)) << 126 << std::setw(3*(6+precVerf)) << " local direction" << " "
118 << G4endl; 127 << G4endl;
119 G4cout << " Mother " << std::setw(12) << 128 G4cout << " Mother " << std::setw(12) << motherSolid->GetName() << " "
120 << std::setw(4+precVerf) << 129 << std::setw(4+precVerf) << localPoint << " "
121 << std::setw(4+precVerf) << 130 << std::setw(4+precVerf) << motherSafety << " "
122 << G4endl; 131 << G4endl;
123 G4cout.precision(oldprec); 132 G4cout.precision(oldprec);
124 } 133 }
125 } 134 }
126 135
127 // ******************************************* 136 // ********************************************************************
128 // AlongComputeStepLog 137 // AlongComputeStepLog
129 // ******************************************* 138 // ********************************************************************
130 // 139 //
131 void 140 void
132 G4NavigationLogger::AlongComputeStepLog(const 141 G4NavigationLogger::AlongComputeStepLog(const G4VSolid* sampleSolid,
133 const 142 const G4ThreeVector& samplePoint,
134 const 143 const G4ThreeVector& sampleDirection,
135 const 144 const G4ThreeVector& localDirection,
136 145 G4double sampleSafety,
137 146 G4double sampleStep) const
138 { 147 {
139 // Check to see that the resulting point is 148 // Check to see that the resulting point is indeed in/on volume.
140 // This check could eventually be made only 149 // This check could eventually be made only for successful candidate.
141 150
142 if ( sampleStep < kInfinity ) 151 if ( sampleStep < kInfinity )
143 { 152 {
144 G4ThreeVector intersectionPoint; 153 G4ThreeVector intersectionPoint;
145 intersectionPoint = samplePoint + sampleSt << 154 intersectionPoint= samplePoint + sampleStep * sampleDirection;
146 EInside insideIntPt = sampleSolid->Inside( << 155 EInside insideIntPt= sampleSolid->Inside(intersectionPoint);
147 G4String fType = fId + "::ComputeStep()"; 156 G4String fType = fId + "::ComputeStep()";
148 157
149 G4String solidResponse = "-kInside-"; 158 G4String solidResponse = "-kInside-";
150 if (insideIntPt == kOutside) 159 if (insideIntPt == kOutside)
151 { solidResponse = "-kOutside-"; } 160 { solidResponse = "-kOutside-"; }
152 else if (insideIntPt == kSurface) 161 else if (insideIntPt == kSurface)
153 { solidResponse = "-kSurface-"; } 162 { solidResponse = "-kSurface-"; }
154 163
155 if ( fVerbose == 1 || fVerbose > 4 ) 164 if ( fVerbose == 1 || fVerbose > 4 )
156 { 165 {
157 G4cout << " Invoked Inside() for soli 166 G4cout << " Invoked Inside() for solid: "
158 << sampleSolid->GetName() 167 << sampleSolid->GetName()
159 << ". Solid replied: " << solidRe 168 << ". Solid replied: " << solidResponse << G4endl
160 << " For point p: " << interse 169 << " For point p: " << intersectionPoint
161 << ", considered as 'intersection 170 << ", considered as 'intersection' point." << G4endl;
162 } 171 }
163 172
164 G4double safetyIn = -1, safetyOut = -1; / << 173 G4double safetyIn= -1, safetyOut= -1; // Set to invalid values
165 G4double newDistIn = -1, newDistOut = -1; << 174 G4double newDistIn= -1, newDistOut= -1;
166 if( insideIntPt != kInside ) 175 if( insideIntPt != kInside )
167 { 176 {
168 safetyIn = sampleSolid->DistanceToIn(int << 177 safetyIn= sampleSolid->DistanceToIn(intersectionPoint);
169 newDistIn = sampleSolid->DistanceToIn(in << 178 newDistIn= sampleSolid->DistanceToIn(intersectionPoint,
170 sa << 179 sampleDirection);
171 } 180 }
172 if( insideIntPt != kOutside ) 181 if( insideIntPt != kOutside )
173 { 182 {
174 safetyOut = sampleSolid->DistanceToOut(i << 183 safetyOut= sampleSolid->DistanceToOut(intersectionPoint);
175 newDistOut = sampleSolid->DistanceToOut( << 184 newDistOut= sampleSolid->DistanceToOut(intersectionPoint,
176 << 185 sampleDirection);
177 } 186 }
178 if( insideIntPt != kSurface ) 187 if( insideIntPt != kSurface )
179 { 188 {
180 std::ostringstream message; 189 std::ostringstream message;
181 message.precision(16); 190 message.precision(16);
182 message << "Conflicting response from So 191 message << "Conflicting response from Solid." << G4endl
183 << " Inaccurate solid D 192 << " Inaccurate solid DistanceToIn"
184 << " for solid " << sampleSolid- 193 << " for solid " << sampleSolid->GetName() << G4endl
185 << " Solid gave Distanc 194 << " Solid gave DistanceToIn = "
186 << sampleStep << " yet returns " 195 << sampleStep << " yet returns " << solidResponse
187 << " for this point !" << G4endl 196 << " for this point !" << G4endl
188 << " Original Point 197 << " Original Point = " << samplePoint << G4endl
189 << " Original Direction 198 << " Original Direction = " << sampleDirection << G4endl
190 << " Intersection Point 199 << " Intersection Point = " << intersectionPoint << G4endl
191 << " Safety values: " 200 << " Safety values: " << G4endl;
192 if ( insideIntPt != kInside ) 201 if ( insideIntPt != kInside )
193 { 202 {
194 message << " DistanceToIn(p) 203 message << " DistanceToIn(p) = " << safetyIn;
195 } 204 }
196 if ( insideIntPt != kOutside ) 205 if ( insideIntPt != kOutside )
197 { 206 {
198 message << " DistanceToOut(p) 207 message << " DistanceToOut(p) = " << safetyOut;
199 } 208 }
200 message << G4endl; 209 message << G4endl;
201 message << " Solid Parameters: " << *sam 210 message << " Solid Parameters: " << *sampleSolid;
202 G4Exception(fType, "GeomNav1001", JustWa 211 G4Exception(fType, "GeomNav1001", JustWarning, message);
203 } 212 }
204 else 213 else
205 { 214 {
206 // If it is on the surface, *ensure* tha 215 // If it is on the surface, *ensure* that either DistanceToIn
207 // or DistanceToOut returns a finite val 216 // or DistanceToOut returns a finite value ( >= Tolerance).
208 // 217 //
209 if( std::max( newDistIn, newDistOut ) <= 218 if( std::max( newDistIn, newDistOut ) <=
210 G4GeometryTolerance::GetInstance()-> 219 G4GeometryTolerance::GetInstance()->GetSurfaceTolerance() )
211 { 220 {
212 std::ostringstream message; 221 std::ostringstream message;
213 message << "Zero from both Solid Dista 222 message << "Zero from both Solid DistanceIn and Out(p,v)." << G4endl
214 << " Identified point for whi 223 << " Identified point for which the solid "
215 << sampleSolid->GetName() << G 224 << sampleSolid->GetName() << G4endl
216 << " has MAJOR problem: " << 225 << " has MAJOR problem: " << G4endl
217 << " --> Both DistanceToIn(p, 226 << " --> Both DistanceToIn(p,v) and DistanceToOut(p,v) "
218 << "return Zero, an equivalent 227 << "return Zero, an equivalent value or negative value."
219 << G4endl 228 << G4endl
220 << " Solid: " << sampleSoli 229 << " Solid: " << sampleSolid << G4endl
221 << " Point p= " << intersec 230 << " Point p= " << intersectionPoint << G4endl
222 << " Direction v= " << samp 231 << " Direction v= " << sampleDirection << G4endl
223 << " DistanceToIn(p,v) 232 << " DistanceToIn(p,v) = " << newDistIn << G4endl
224 << " DistanceToOut(p,v,..) 233 << " DistanceToOut(p,v,..) = " << newDistOut << G4endl
225 << " Safety values: " << G4 234 << " Safety values: " << G4endl
226 << " DistanceToIn(p) = " 235 << " DistanceToIn(p) = " << safetyIn << G4endl
227 << " DistanceToOut(p) = " 236 << " DistanceToOut(p) = " << safetyOut;
228 G4Exception(fType, "GeomNav0003", Fata 237 G4Exception(fType, "GeomNav0003", FatalException, message);
229 } 238 }
230 } 239 }
231 240
232 // Verification / verbosity 241 // Verification / verbosity
233 // 242 //
234 if ( fVerbose > 1 ) 243 if ( fVerbose > 1 )
235 { 244 {
236 static const G4int precVerf= 20; // Pre 245 static const G4int precVerf= 20; // Precision
237 G4long oldprec = G4cout.precision(precVe << 246 G4int oldprec = G4cout.precision(precVerf);
238 G4cout << "Daughter " 247 G4cout << "Daughter "
239 << std::setw(12) << sampl 248 << std::setw(12) << sampleSolid->GetName() << " "
240 << std::setw(4+precVerf) << sampl 249 << std::setw(4+precVerf) << samplePoint << " "
241 << std::setw(4+precVerf) << sampl 250 << std::setw(4+precVerf) << sampleSafety << " "
242 << std::setw(4+precVerf) << sampl 251 << std::setw(4+precVerf) << sampleStep << " "
243 << std::setw(16) << "dist 252 << std::setw(16) << "distanceToIn" << " "
244 << std::setw(4+precVerf) << local 253 << std::setw(4+precVerf) << localDirection << " "
245 << G4endl; 254 << G4endl;
246 G4cout.precision(oldprec); 255 G4cout.precision(oldprec);
247 } 256 }
248 } 257 }
249 } 258 }
250 259
251 // ******************************************* 260 // ********************************************************************
252 // CheckDaughterEntryPoint 261 // CheckDaughterEntryPoint
253 // ******************************************* 262 // ********************************************************************
254 // 263 //
255 void 264 void
256 G4NavigationLogger::CheckDaughterEntryPoint(co 265 G4NavigationLogger::CheckDaughterEntryPoint(const G4VSolid* sampleSolid,
257 co 266 const G4ThreeVector& samplePoint,
258 co 267 const G4ThreeVector& sampleDirection,
259 co 268 const G4VSolid* motherSolid,
260 co 269 const G4ThreeVector& localPoint,
261 co 270 const G4ThreeVector& localDirection,
262 271 G4double motherStep,
263 272 G4double sampleStep) const
264 { 273 {
265 const G4double kCarTolerance = motherSolid-> 274 const G4double kCarTolerance = motherSolid->GetTolerance();
266 275
267 // Double check the expected condition of be 276 // Double check the expected condition of being called
268 // 277 //
269 G4bool SuspiciousDaughterDist = ( sampleStep 278 G4bool SuspiciousDaughterDist = ( sampleStep >= motherStep )
270 && ( sampleStep 279 && ( sampleStep < kInfinity );
271 280
272 if( sampleStep >= kInfinity ) 281 if( sampleStep >= kInfinity )
273 { 282 {
274 G4ExceptionDescription msg; 283 G4ExceptionDescription msg;
275 msg.precision(12); 284 msg.precision(12);
276 msg << " WARNING - Called with 'infinite' 285 msg << " WARNING - Called with 'infinite' step. " << G4endl;
277 msg << " Checks have no meaning if daug 286 msg << " Checks have no meaning if daughter step is infinite." << G4endl;
278 msg << " kInfinity = " << kInfinity / 287 msg << " kInfinity = " << kInfinity / millimeter << G4endl;
279 msg << " sampleStep = " << sampleStep / 288 msg << " sampleStep = " << sampleStep / millimeter << G4endl;
280 msg << " sampleStep < kInfinity " << (s 289 msg << " sampleStep < kInfinity " << (sampleStep<kInfinity) << G4endl;
281 msg << " kInfinity - sampleStep " << (k 290 msg << " kInfinity - sampleStep " << (kInfinity-sampleStep) / millimeter << G4endl;
282 msg << " Returning immediately."; 291 msg << " Returning immediately.";
283 G4Exception("G4NavigationLogger::CheckDaug 292 G4Exception("G4NavigationLogger::CheckDaughterEntryPoint()",
284 "GeomNav0003", JustWarning, ms 293 "GeomNav0003", JustWarning, msg);
285 return; 294 return;
286 } 295 }
287 296
288 // The intersection point with the daughter 297 // The intersection point with the daughter is after the exit point
289 // from the mother volume !! 298 // from the mother volume !!
290 // This is legal / allowed to occur only if 299 // This is legal / allowed to occur only if the mother is concave
291 // ***************************************** 300 // ****************************************************************
292 // If mother is convex the daughter volume m 301 // If mother is convex the daughter volume must be encountered
293 // before the exit from the current volume! 302 // before the exit from the current volume!
294 303
295 // Check #1) whether the track will re-enter 304 // Check #1) whether the track will re-enter the current mother
296 // in the extension past its curre 305 // in the extension past its current exit point
297 G4ThreeVector localExitMotherPos = localPoin << 306 G4ThreeVector localExitMotherPos= localPoint+motherStep*localDirection;
298 G4double distExitToReEntry = motherSolid->Di << 307 G4double distExitToReEntry= motherSolid->DistanceToIn(localExitMotherPos,
299 << 308 localDirection);
300 309
301 // Check #2) whether the 'entry' point in th 310 // Check #2) whether the 'entry' point in the daughter is inside the mother
302 // 311 //
303 G4ThreeVector localEntryInDaughter = localPo 312 G4ThreeVector localEntryInDaughter = localPoint+sampleStep*localDirection;
304 EInside insideMother = motherSolid->Inside( << 313 EInside insideMother= motherSolid->Inside( localEntryInDaughter );
305 314
306 G4String solidResponse = "-kInside-"; 315 G4String solidResponse = "-kInside-";
307 if (insideMother == kOutside) { solidR 316 if (insideMother == kOutside) { solidResponse = "-kOutside-"; }
308 else if (insideMother == kSurface) { solidR 317 else if (insideMother == kSurface) { solidResponse = "-kSurface-"; }
309 318
310 G4double distToReEntry = distExitToReE 319 G4double distToReEntry = distExitToReEntry + motherStep;
311 G4ThreeVector localReEntryPoint = localPoin 320 G4ThreeVector localReEntryPoint = localPoint+distToReEntry*localDirection;
312 321
313 // Clear error -- Daughter entry point is b 322 // Clear error -- Daughter entry point is bad
314 constexpr G4double eps= 1.0e-10; <<
315 G4bool DaughterEntryIsOutside = SuspiciousDa 323 G4bool DaughterEntryIsOutside = SuspiciousDaughterDist
316 && ( (sampleStep * (1.0+eps) < distToReEn << 324 && ( (sampleStep < distToReEntry) || (insideMother == kOutside ) );
317 G4bool EntryIsMotherExit = std::fabs(sampleS 325 G4bool EntryIsMotherExit = std::fabs(sampleStep-motherStep) < kCarTolerance;
318 326
319 // Check for more subtle error - is exit poi 327 // Check for more subtle error - is exit point of daughter correct ?
320 G4ThreeVector sampleEntryPoint = samplePoint << 328 G4ThreeVector sampleEntryPoint= samplePoint+sampleStep*sampleDirection;
321 G4double sampleCrossingDist = sampleSolid->D << 329 G4double sampleCrossingDist= sampleSolid->DistanceToOut( sampleEntryPoint,
322 << 330 sampleDirection );
323 G4double sampleExitDist = sampleStep+sa 331 G4double sampleExitDist = sampleStep+sampleCrossingDist;
324 G4ThreeVector sampleExitPoint = samplePoint+ << 332 G4ThreeVector sampleExitPoint= samplePoint+sampleExitDist*sampleDirection;
325 333
326 G4bool TransitProblem = ( (sampleStep < moth 334 G4bool TransitProblem = ( (sampleStep < motherStep)
327 && (sampleExitDist > 335 && (sampleExitDist > motherStep + kCarTolerance) )
328 || ( EntryIsMotherExit && (sampleCr 336 || ( EntryIsMotherExit && (sampleCrossingDist > kCarTolerance) );
329 337
330 if( DaughterEntryIsOutside 338 if( DaughterEntryIsOutside
331 || TransitProblem 339 || TransitProblem
332 || (SuspiciousDaughterDist && (fVerbose 340 || (SuspiciousDaughterDist && (fVerbose > 3) ) )
333 { 341 {
334 G4ExceptionDescription msg; 342 G4ExceptionDescription msg;
335 msg.precision(16); 343 msg.precision(16);
336 344
337 if( DaughterEntryIsOutside ) 345 if( DaughterEntryIsOutside )
338 { 346 {
339 msg << "WARNING> Intersection distance t 347 msg << "WARNING> Intersection distance to Daughter volume is further"
340 << " than the distance to boundar 348 << " than the distance to boundary." << G4endl
341 << " It appears that part of the da 349 << " It appears that part of the daughter volume is *outside*"
342 << " this mother. " << G4endl; 350 << " this mother. " << G4endl;
343 msg << " One of the following checks si 351 msg << " One of the following checks signaled a problem:" << G4endl
344 << " -sampleStep (dist to daugh) < 352 << " -sampleStep (dist to daugh) < mother-exit dist + distance "
345 << "to ReEntry point for mother 353 << "to ReEntry point for mother " << G4endl
346 << " -position of daughter intersec 354 << " -position of daughter intersection is outside mother volume."
347 << G4endl; 355 << G4endl;
348 } 356 }
349 else if( TransitProblem ) 357 else if( TransitProblem )
350 { 358 {
351 G4double protrusion = sampleExitDist - m 359 G4double protrusion = sampleExitDist - motherStep;
352 360
353 msg << "WARNING> Daughter volume extend 361 msg << "WARNING> Daughter volume extends beyond mother boundary. "
354 << G4endl; 362 << G4endl;
355 if ( ( sampleStep < motherStep ) 363 if ( ( sampleStep < motherStep )
356 && (sampleExitDist > motherStep + kCar 364 && (sampleExitDist > motherStep + kCarTolerance ) )
357 { 365 {
358 // 1st Issue with Daughter 366 // 1st Issue with Daughter
359 msg << " Crossing distance in t 367 msg << " Crossing distance in the daughter causes is to extend"
360 << " beyond the mother exit. " << 368 << " beyond the mother exit. " << G4endl;
361 msg << " Length protruding = " 369 msg << " Length protruding = " << protrusion << G4endl;
362 } 370 }
363 if( EntryIsMotherExit ) 371 if( EntryIsMotherExit )
364 { 372 {
365 // 1st Issue with Daughter 373 // 1st Issue with Daughter
366 msg << " Intersection distance 374 msg << " Intersection distance to Daughter is within "
367 << " tolerance of the distance" << 375 << " tolerance of the distance" << G4endl;
368 msg << " to the mother boundary 376 msg << " to the mother boundary * and * " << G4endl;
369 msg << " the crossing distance 377 msg << " the crossing distance in the daughter is > tolerance."
370 << G4endl; 378 << G4endl;
371 } 379 }
372 } 380 }
373 else 381 else
374 { 382 {
375 msg << "NearMiss> Intersection to Daught 383 msg << "NearMiss> Intersection to Daughter volume is in extension past the"
376 << " current exit point of the mot 384 << " current exit point of the mother volume." << G4endl;
377 msg << " This is not an error - << 385 msg << " This is not an error - just an unusual occurence,"
378 << " possible in the case of conca 386 << " possible in the case of concave volume. " << G4endl;
379 } 387 }
380 msg << "---- Information about intersectio 388 msg << "---- Information about intersection with daughter, mother: "
381 << G4endl; 389 << G4endl;
382 msg << " sampleStep (daughter) = " << s 390 msg << " sampleStep (daughter) = " << sampleStep << G4endl
383 << " motherStep = " << m 391 << " motherStep = " << motherStep << G4endl
384 << " distToRentry(mother) = " << d 392 << " distToRentry(mother) = " << distToReEntry << G4endl
385 << " Inside(entry pnt daug): " << s 393 << " Inside(entry pnt daug): " << solidResponse << G4endl
386 << " dist across daughter = " << s 394 << " dist across daughter = " << sampleCrossingDist << G4endl;
387 msg << " Mother Name (Solid) : " << mother 395 msg << " Mother Name (Solid) : " << motherSolid->GetName() << G4endl
388 << " In local (mother) coordinates: " 396 << " In local (mother) coordinates: " << G4endl
389 << " Starting Point = " << l 397 << " Starting Point = " << localPoint << G4endl
390 << " Direction = " << l 398 << " Direction = " << localDirection << G4endl
391 << " Exit Point (mother)= " << l 399 << " Exit Point (mother)= " << localExitMotherPos << G4endl
392 << " Entry Point (daughter)= " << l 400 << " Entry Point (daughter)= " << localPoint+sampleStep*localDirection
393 << G4endl; 401 << G4endl;
394 if( distToReEntry < kInfinity ) 402 if( distToReEntry < kInfinity )
395 { 403 {
396 msg << " ReEntry Point (mother)= " << 404 msg << " ReEntry Point (mother)= " << localReEntryPoint << G4endl;
397 } 405 }
398 else 406 else
399 { 407 {
400 msg << " No ReEntry - track does not 408 msg << " No ReEntry - track does not encounter mother volume again! "
401 << G4endl; 409 << G4endl;
402 } 410 }
403 msg << " Daughter Name (Solid): " << sampl 411 msg << " Daughter Name (Solid): " << sampleSolid->GetName() << G4endl
404 << " In daughter coordinates: " << G4e 412 << " In daughter coordinates: " << G4endl
405 << " Starting Point = " << s 413 << " Starting Point = " << samplePoint << G4endl
406 << " Direction = " << s 414 << " Direction = " << sampleDirection << G4endl
407 << " Entry Point (daughter)= " << s 415 << " Entry Point (daughter)= " << sampleEntryPoint
408 << G4endl; 416 << G4endl;
409 msg << " Description of mother solid: " < 417 msg << " Description of mother solid: " << G4endl
410 << *motherSolid << G4endl 418 << *motherSolid << G4endl
411 << " Description of daughter solid: " 419 << " Description of daughter solid: " << G4endl
412 << *sampleSolid << G4endl; 420 << *sampleSolid << G4endl;
413 G4String fType = fId + "::ComputeStep()"; 421 G4String fType = fId + "::ComputeStep()";
414 422
415 if( DaughterEntryIsOutside || TransitProbl 423 if( DaughterEntryIsOutside || TransitProblem )
416 { 424 {
417 G4Exception(fType, "GeomNav0003", JustWa 425 G4Exception(fType, "GeomNav0003", JustWarning, msg);
418 } 426 }
419 else 427 else
420 { 428 {
421 G4cout << fType 429 G4cout << fType
422 << " -- Checked distance of Entry 430 << " -- Checked distance of Entry to daughter vs exit of mother"
423 << G4endl; 431 << G4endl;
424 G4cout << msg.str(); 432 G4cout << msg.str();
425 G4cout << G4endl; 433 G4cout << G4endl;
426 } 434 }
427 } 435 }
428 } 436 }
429 437
430 // ******************************************* 438 // ********************************************************************
431 // PostComputeStepLog 439 // PostComputeStepLog
432 // ******************************************* 440 // ********************************************************************
433 // 441 //
434 void 442 void
435 G4NavigationLogger::PostComputeStepLog(const G 443 G4NavigationLogger::PostComputeStepLog(const G4VSolid* motherSolid,
436 const G 444 const G4ThreeVector& localPoint,
437 const G 445 const G4ThreeVector& localDirection,
438 G 446 G4double motherStep,
439 G 447 G4double motherSafety) const
440 { 448 {
441 if ( fVerbose == 1 || fVerbose > 4 ) 449 if ( fVerbose == 1 || fVerbose > 4 )
442 { 450 {
443 G4cout << " Mother " 451 G4cout << " Mother "
444 << std::setw(15) << motherSafety << 452 << std::setw(15) << motherSafety << " "
445 << std::setw(15) << motherStep << 453 << std::setw(15) << motherStep << " " << localPoint << " - "
446 << motherSolid->GetEntityType() << 454 << motherSolid->GetEntityType() << ": " << motherSolid->GetName()
447 << G4endl; 455 << G4endl;
448 } 456 }
449 if( ( motherStep < 0.0 ) || ( motherStep >= 457 if( ( motherStep < 0.0 ) || ( motherStep >= kInfinity) )
450 { 458 {
451 G4String fType = fId + "::ComputeStep()"; 459 G4String fType = fId + "::ComputeStep()";
452 G4long oldPrOut = G4cout.precision(16); << 460 G4int oldPrOut= G4cout.precision(16);
453 G4long oldPrErr = G4cerr.precision(16); << 461 G4int oldPrErr= G4cerr.precision(16);
454 std::ostringstream message; 462 std::ostringstream message;
455 message << "Current point is outside the c 463 message << "Current point is outside the current solid !" << G4endl
456 << " Problem in Navigation" 464 << " Problem in Navigation" << G4endl
457 << " Point (local coordinat 465 << " Point (local coordinates): "
458 << localPoint << G4endl 466 << localPoint << G4endl
459 << " Local Direction: " << 467 << " Local Direction: " << localDirection << G4endl
460 << " Solid: " << motherSoli 468 << " Solid: " << motherSolid->GetName();
461 motherSolid->DumpInfo(); 469 motherSolid->DumpInfo();
462 G4Exception(fType, "GeomNav0003", FatalExc 470 G4Exception(fType, "GeomNav0003", FatalException, message);
463 G4cout.precision(oldPrOut); 471 G4cout.precision(oldPrOut);
464 G4cerr.precision(oldPrErr); 472 G4cerr.precision(oldPrErr);
465 } 473 }
466 if ( fVerbose > 1 ) 474 if ( fVerbose > 1 )
467 { 475 {
468 static const G4int precVerf = 20; // Prec << 476 static const G4int precVerf= 20; // Precision
469 G4long oldprec = G4cout.precision(precVerf << 477 G4int oldprec = G4cout.precision(precVerf);
470 G4cout << " Mother " << std::setw(12) << 478 G4cout << " Mother " << std::setw(12) << motherSolid->GetName() << " "
471 << std::setw(4+precVerf) << l 479 << std::setw(4+precVerf) << localPoint << " "
472 << std::setw(4+precVerf) << m 480 << std::setw(4+precVerf) << motherSafety << " "
473 << std::setw(4+precVerf) << m 481 << std::setw(4+precVerf) << motherStep << " "
474 << std::setw(16) << " 482 << std::setw(16) << "distanceToOut" << " "
475 << std::setw(4+precVerf) << l 483 << std::setw(4+precVerf) << localDirection << " "
476 << G4endl; 484 << G4endl;
477 G4cout.precision(oldprec); 485 G4cout.precision(oldprec);
478 } 486 }
479 } 487 }
480 488
481 // ******************************************* 489 // ********************************************************************
482 // ComputeSafetyLog 490 // ComputeSafetyLog
483 // ******************************************* 491 // ********************************************************************
484 // 492 //
485 void 493 void
486 G4NavigationLogger::ComputeSafetyLog(const G4V 494 G4NavigationLogger::ComputeSafetyLog(const G4VSolid* solid,
487 const G4T 495 const G4ThreeVector& point,
488 G4d 496 G4double safety,
489 G4b 497 G4bool isMotherVolume,
490 G4i 498 G4int banner) const
491 { 499 {
492 if( banner < 0 ) 500 if( banner < 0 )
493 { 501 {
494 banner = static_cast<G4int>(isMotherVolume << 502 banner = isMotherVolume;
495 } 503 }
496 if( fVerbose >= 1 ) 504 if( fVerbose >= 1 )
497 { 505 {
498 G4String volumeType = isMotherVolume ? " M 506 G4String volumeType = isMotherVolume ? " Mother " : "Daughter";
499 if (banner != 0) << 507 if (banner)
500 { 508 {
501 G4cout << "************** " << fId << ": 509 G4cout << "************** " << fId << "::ComputeSafety() ****************"
502 << G4endl; 510 << G4endl;
503 G4cout << " VolType " 511 G4cout << " VolType "
504 << std::setw(15) << "Safety/mm" < 512 << std::setw(15) << "Safety/mm" << " "
505 << std::setw(52) << "Position (lo 513 << std::setw(52) << "Position (local coordinates)"
506 << " - Solid" << G4endl; 514 << " - Solid" << G4endl;
507 } 515 }
508 G4cout << volumeType 516 G4cout << volumeType
509 << std::setw(15) << safety << " " < 517 << std::setw(15) << safety << " " << point << " - "
510 << solid->GetEntityType() << ": " < 518 << solid->GetEntityType() << ": " << solid->GetName() << G4endl;
511 } 519 }
512 } 520 }
513 521
514 // ******************************************* 522 // ********************************************************************
515 // PrintDaughterLog 523 // PrintDaughterLog
516 // ******************************************* 524 // ********************************************************************
517 // 525 //
518 void 526 void
519 G4NavigationLogger::PrintDaughterLog (const G4 527 G4NavigationLogger::PrintDaughterLog (const G4VSolid* sampleSolid,
520 const G4 528 const G4ThreeVector& samplePoint,
521 G4 529 G4double sampleSafety,
522 G4 530 G4bool withStep,
523 const G4 531 const G4ThreeVector& sampleDirection, G4double sampleStep ) const
524 { 532 {
525 if ( fVerbose >= 1 ) 533 if ( fVerbose >= 1 )
526 { 534 {
527 G4long oldPrec = G4cout.precision(8); << 535 G4int oldPrec= G4cout.precision(8);
528 G4cout << "Daughter " 536 G4cout << "Daughter "
529 << std::setw(15) << sampleSafety << 537 << std::setw(15) << sampleSafety << " ";
530 if (withStep) // (sampleStep != -1.0 ) 538 if (withStep) // (sampleStep != -1.0 )
531 { 539 {
532 G4cout << std::setw(15) << sampleStep << 540 G4cout << std::setw(15) << sampleStep << " ";
533 } 541 }
534 else 542 else
535 { 543 {
536 G4cout << std::setw(15) << "Not-Availabl 544 G4cout << std::setw(15) << "Not-Available" << " ";
537 } 545 }
538 G4cout << samplePoint << " - " 546 G4cout << samplePoint << " - "
539 << sampleSolid->GetEntityType() << 547 << sampleSolid->GetEntityType() << ": " << sampleSolid->GetName();
540 if( withStep ) 548 if( withStep )
541 { 549 {
542 G4cout << " dir= " << sampleDirection; 550 G4cout << " dir= " << sampleDirection;
543 } 551 }
544 G4cout << G4endl; 552 G4cout << G4endl;
545 G4cout.precision(oldPrec); 553 G4cout.precision(oldPrec);
546 } 554 }
547 } 555 }
548 556
549 // ******************************************* 557 // ********************************************************************
550 // CheckAndReportBadNormal 558 // CheckAndReportBadNormal
551 // ******************************************* 559 // ********************************************************************
552 // 560 //
553 G4bool 561 G4bool
554 G4NavigationLogger:: 562 G4NavigationLogger::
555 CheckAndReportBadNormal(const G4ThreeVector& u 563 CheckAndReportBadNormal(const G4ThreeVector& unitNormal,
556 const G4ThreeVector& l 564 const G4ThreeVector& localPoint,
557 const G4ThreeVector& l 565 const G4ThreeVector& localDirection,
558 G4double 566 G4double step,
559 const G4VSolid* s 567 const G4VSolid* solid,
560 const char* m 568 const char* msg ) const
561 { 569 {
562 G4double normMag2 = unitNormal.mag2(); << 570 G4double normMag2 = unitNormal.mag2();
563 G4bool badLength = ( std::fabs ( normMag2 - 571 G4bool badLength = ( std::fabs ( normMag2 - 1.0 ) > CLHEP::perMillion );
564 572
565 if( badLength ) 573 if( badLength )
566 { 574 {
567 G4double normMag = std::sqrt(normMag2); << 575 G4double normMag= std::sqrt(normMag2);
568 G4ExceptionDescription message; 576 G4ExceptionDescription message;
569 message.precision(10); 577 message.precision(10);
570 message << "============================== 578 message << "============================================================"
571 << G4endl; 579 << G4endl;
572 message << " WARNING> Normal is not a uni 580 message << " WARNING> Normal is not a unit vector. "
573 << " - but |normal| = " << nor 581 << " - but |normal| = " << normMag
574 << " - and |normal|^2 = " << 582 << " - and |normal|^2 = " << normMag2 << G4endl
575 << " which differ from 1.0 by: 583 << " which differ from 1.0 by: " << G4endl
576 << " |normal|-1 = " << norm 584 << " |normal|-1 = " << normMag - 1.0
577 << " and |normal|^2 - 1 = " << 585 << " and |normal|^2 - 1 = " << normMag2 - 1.0 << G4endl
578 << " n = " << unitNormal << G4en 586 << " n = " << unitNormal << G4endl;
579 message << " Info string: " << msg << G4en 587 message << " Info string: " << msg << G4endl;
580 message << "============================== 588 message << "============================================================"
581 << G4endl; 589 << G4endl;
582 590
583 message.precision(16); 591 message.precision(16);
584 592
585 message << " Information on call to Distan 593 message << " Information on call to DistanceToOut: " << G4endl;
586 message << " Position = " << localPoint 594 message << " Position = " << localPoint << G4endl
587 << " Direction = " << localDirec 595 << " Direction = " << localDirection << G4endl;
588 message << " Obtained> distance = " 596 message << " Obtained> distance = " << step << G4endl;
589 message << " > Exit position = " 597 message << " > Exit position = " << localPoint+step*localDirection
590 << G4endl; 598 << G4endl;
591 message << " Parameters of solid: " << 599 message << " Parameters of solid: " << G4endl;
592 message << *solid; 600 message << *solid;
593 message << "============================== 601 message << "============================================================";
594 602
595 G4String fMethod = fId + "::ComputeStep()" 603 G4String fMethod = fId + "::ComputeStep()";
596 G4Exception( fMethod, "GeomNav0003", JustW 604 G4Exception( fMethod, "GeomNav0003", JustWarning, message);
597 } 605 }
598 return badLength; 606 return badLength;
599 } 607 }
600 608
601 // ******************************************* 609 // ********************************************************************
602 // CheckAndReportBadNormal - due to Rotation M 610 // CheckAndReportBadNormal - due to Rotation Matrix
603 // ******************************************* 611 // ********************************************************************
604 // 612 //
605 G4bool 613 G4bool
606 G4NavigationLogger:: 614 G4NavigationLogger::
607 CheckAndReportBadNormal(const G4ThreeVector& r 615 CheckAndReportBadNormal(const G4ThreeVector& rotatedNormal,
608 const G4ThreeVector& o 616 const G4ThreeVector& originalNormal,
609 const G4RotationMatrix 617 const G4RotationMatrix& rotationM,
610 const char* m 618 const char* msg ) const
611 { 619 {
612 G4double normMag2 = rotatedNormal.mag2(); 620 G4double normMag2 = rotatedNormal.mag2();
613 G4bool badLength = ( std::fabs ( normMag2 - 621 G4bool badLength = ( std::fabs ( normMag2 - 1.0 ) > CLHEP::perMillion );
614 622
615 if( badLength ) 623 if( badLength )
616 { 624 {
617 G4double normMag = std::sqrt(normMag2); << 625 G4double normMag= std::sqrt(normMag2);
618 G4ExceptionDescription message; 626 G4ExceptionDescription message;
619 message.precision(10); 627 message.precision(10);
620 message << "============================== 628 message << "============================================================"
621 << G4endl; 629 << G4endl;
622 message << " WARNING> Rotated n(ormal) is 630 message << " WARNING> Rotated n(ormal) is not a unit vector. " << G4endl
623 << " |normal| = " << normMa 631 << " |normal| = " << normMag
624 << " and |normal|^2 = " << 632 << " and |normal|^2 = " << normMag2 << G4endl
625 << " Diff from 1.0: " << G4endl 633 << " Diff from 1.0: " << G4endl
626 << " |normal|-1 = " << normMag 634 << " |normal|-1 = " << normMag - 1.0
627 << " and |normal|^2 - 1 = " << n 635 << " and |normal|^2 - 1 = " << normMag2 - 1.0 << G4endl;
628 message << " Rotated n = (" << rotatedN 636 message << " Rotated n = (" << rotatedNormal.x() << "," << rotatedNormal.y() << ","
629 << rotatedNormal.z() << ")" << G4e 637 << rotatedNormal.z() << ")" << G4endl;
630 message << " Original n = (" << original 638 message << " Original n = (" << originalNormal.x() << "," << originalNormal.y() << ","
631 << originalNormal.z() << ")" << G4 639 << originalNormal.z() << ")" << G4endl;
632 message << " Info string: " << msg << G4en 640 message << " Info string: " << msg << G4endl;
633 message << "============================== 641 message << "============================================================"
634 << G4endl; 642 << G4endl;
635 643
636 message.precision(16); 644 message.precision(16);
637 645
638 message << " Information on RotationMatrix 646 message << " Information on RotationMatrix : " << G4endl;
639 message << " Original: " << G4endl; 647 message << " Original: " << G4endl;
640 message << rotationM << G4endl; 648 message << rotationM << G4endl;
641 message << " Inverse (used in transformati 649 message << " Inverse (used in transformation): " << G4endl;
642 message << rotationM.inverse() << G4endl; 650 message << rotationM.inverse() << G4endl;
643 message << "============================== 651 message << "============================================================";
644 652
645 G4String fMethod = fId + "::ComputeStep()" 653 G4String fMethod = fId + "::ComputeStep()";
646 G4Exception( fMethod, "GeomNav0003", JustW 654 G4Exception( fMethod, "GeomNav0003", JustWarning, message);
647 } 655 }
648 return badLength; 656 return badLength;
649 } 657 }
650 658
651 // ******************************************* 659 // ********************************************************************
652 // ReportOutsideMother 660 // ReportOutsideMother
653 // ******************************************* 661 // ********************************************************************
654 // 662 //
655 // Report Exception if point is outside mother 663 // Report Exception if point is outside mother.
656 // Fatal exception will be used if either 'che 664 // Fatal exception will be used if either 'checkMode error is > triggerDist
657 // 665 //
658 void 666 void
659 G4NavigationLogger::ReportOutsideMother(const 667 G4NavigationLogger::ReportOutsideMother(const G4ThreeVector& localPoint,
660 const 668 const G4ThreeVector& localDirection,
661 const 669 const G4VPhysicalVolume* physical,
662 670 G4double triggerDist) const
663 { 671 {
664 const G4LogicalVolume* logicalVol = physical << 672 const G4LogicalVolume* logicalVol = physical
665 ? physical << 673 ? physical->GetLogicalVolume() : 0;
666 const G4VSolid* solid = logicalVol != nullpt << 674 const G4VSolid* solid = logicalVol
667 ? logicalVol->GetSolid << 675 ? logicalVol->GetSolid() : 0;
668 676
669 G4String fMethod = fId + "::ComputeStep()"; 677 G4String fMethod = fId + "::ComputeStep()";
670 678
671 if( solid == nullptr ) << 679 if( solid == 0 )
672 { 680 {
673 G4Exception(fMethod, "GeomNav0003", FatalE 681 G4Exception(fMethod, "GeomNav0003", FatalException,
674 "Erroneous call to ReportOutsi 682 "Erroneous call to ReportOutsideMother: no Solid is available");
675 return; 683 return;
676 } 684 }
677 const G4double kCarTolerance = solid->GetTol 685 const G4double kCarTolerance = solid->GetTolerance();
678 686
679 // Double check reply - it should be kInfini 687 // Double check reply - it should be kInfinity
680 const G4double distToOut = solid->DistanceTo << 688 const G4double distToOut = solid->DistanceToOut(localPoint, localDirection);
>> 689 // const G4double distToOutNeg = solid->DistanceToOut(localPoint, -localDirection);
>> 690
681 const EInside inSolid = solid->Inside(loc 691 const EInside inSolid = solid->Inside(localPoint);
682 const G4double safetyToIn = solid->Distance 692 const G4double safetyToIn = solid->DistanceToIn(localPoint);
683 const G4double safetyToOut = solid->Distance 693 const G4double safetyToOut = solid->DistanceToOut(localPoint);
684 // const G4double distToInPos = << 694
685 // solid->DistanceToIn(localP << 695 const G4double distToInPos = solid->DistanceToIn(localPoint, localDirection);
>> 696 // const G4double distToInNeg = solid->DistanceToIn(localPoint, -localDirection);
686 697
687 // 1. Check consistency between Safety obtai 698 // 1. Check consistency between Safety obtained and report from distance
688 // We must ensure that (mother)Safety <= 699 // We must ensure that (mother)Safety <= 0.0
689 // in the case that the point is outsi 700 // in the case that the point is outside!
690 // [ If this is not the case, this proble 701 // [ If this is not the case, this problem can easily go undetected,
691 // except in Check mode ! ] 702 // except in Check mode ! ]
692 if( safetyToOut > kCarTolerance 703 if( safetyToOut > kCarTolerance
693 && ( distToOut < 0.0 || distToOut >= kIn 704 && ( distToOut < 0.0 || distToOut >= kInfinity ) )
694 { 705 {
695 G4ExceptionDescription msg1; 706 G4ExceptionDescription msg1;
696 // fNavClerk->ReportBadSafety(localPoint, 707 // fNavClerk->ReportBadSafety(localPoint, localDirection,
697 // motherPhysical, mo 708 // motherPhysical, motherSafety, motherStep);
698 msg1 << " Dangerous inconsistency in resp 709 msg1 << " Dangerous inconsistency in response of solid." << G4endl
699 << " Solid type: " << solid->GetE 710 << " Solid type: " << solid->GetEntityType()
700 << " Name= " << solid->GetName() 711 << " Name= " << solid->GetName() << G4endl;
701 msg1 << " Mother volume gives safety > 0 712 msg1 << " Mother volume gives safety > 0 despite being called for *Outside* point "
702 << G4endl 713 << G4endl
703 << " Location = " << localPoint 714 << " Location = " << localPoint << G4endl
704 << " Direction= " << localDirectio 715 << " Direction= " << localDirection << G4endl
705 << " - Safety (Isotropic d) = " << 716 << " - Safety (Isotropic d) = " << safetyToOut << G4endl
706 << " - Intersection Distance= " << 717 << " - Intersection Distance= " << distToOut << G4endl
707 << G4endl; 718 << G4endl;
708 G4Exception( fMethod, "GeomNav0123", Just 719 G4Exception( fMethod, "GeomNav0123", JustWarning, msg1);
709 } 720 }
710 721
711 // 2. Inconsistency - Too many distances are 722 // 2. Inconsistency - Too many distances are zero (or will be rounded to zero)
712 723
713 // if( std::fabs(distToOut) < kCarTolerance << 724 G4ExceptionDescription msg2;
714 // && std::fabs(distToInPos) < kCarTolerance << 725 if( std::fabs(distToOut) < kCarTolerance && std::fabs(distToInPos) < kCarTolerance )
715 // { << 726 {
716 // If both distanceToIn and distanceToOut << 727 // If both distanceToIn and distanceToOut (p,v) are zero for one direction,
717 // one direction, the particle could get << 728 // the particle could get stuck!
718 // } << 729 }
719 730
720 G4ExceptionDescription msg; 731 G4ExceptionDescription msg;
721 msg.precision(10); 732 msg.precision(10);
>> 733 // G4bool reportIssue= true;
722 734
723 if( std::fabs(distToOut) < kCarTolerance ) 735 if( std::fabs(distToOut) < kCarTolerance )
724 { 736 {
725 // 3. Soft error - safety is not rounded t 737 // 3. Soft error - safety is not rounded to zero
726 // Report nothing - except in 'loud' mo 738 // Report nothing - except in 'loud' mode
727 if( fReportSoftWarnings ) 739 if( fReportSoftWarnings )
728 { 740 {
729 msg << " Warning> DistanceToOut(p,v): " << 741 // reportIssue= true;
730 << "Distance from surface is not rou << 742 msg << " Warning> DistanceToOut(p,v): Distance from surface is not rounded to zero" << G4endl;
731 } << 743 } else {
732 else << 744 // reportIssue= false;
733 { <<
734 return; 745 return;
735 } 746 }
736 } 747 }
737 else 748 else
738 { 749 {
739 // 4. General message - complain that the 750 // 4. General message - complain that the point is outside!
740 // and provide all information about t 751 // and provide all information about the current location,
741 // direction and the answers of the so 752 // direction and the answers of the solid
742 msg << "================================== << 753 msg << "============================================================" << G4endl;
743 << G4endl; << 754 msg << " WARNING> Current Point appears to be Outside mother volume !! " << G4endl;
744 msg << " WARNING> Current Point appears t << 755 msg << " Response of DistanceToOut was negative or kInfinity when called in "
745 << G4endl; << 756 << fMethod << G4endl;
746 msg << " Response of DistanceToOut was n <<
747 << " when called in " << fMethod << G4 <<
748 } 757 }
749 758
750 // Generate and 'print'/stream all the infor 759 // Generate and 'print'/stream all the information needed
751 this->ReportVolumeAndIntersection(msg, local << 760 this->ReportVolumeAndIntersection( msg, localPoint, localDirection, physical );
752 761
753 // Default for distance of 'major' error 762 // Default for distance of 'major' error
754 if( triggerDist <= 0.0 ) << 763 if( triggerDist <= 0.0 ) {
755 { << 764 // triggerDist = 1.e+6 * kCarTolerance; // Well beyond tolerance
756 triggerDist = std::max ( 1.0e+6 * kCarTole 765 triggerDist = std::max ( 1.0e+6 * kCarTolerance, // Well beyond tolerance
757 fMinTriggerDistan 766 fMinTriggerDistance );
758 } 767 }
759 768
760 G4bool majorError = inSolid == kOutside << 769 G4bool majorError = inSolid==kOutside
761 ? ( safetyToIn > triggerDi 770 ? ( safetyToIn > triggerDist )
762 : ( safetyToOut > triggerD 771 : ( safetyToOut > triggerDist );
763 772
764 G4ExceptionSeverity exceptionType = JustWarn 773 G4ExceptionSeverity exceptionType = JustWarning;
765 if ( majorError ) 774 if ( majorError )
766 { 775 {
767 exceptionType = FatalException; 776 exceptionType = FatalException;
768 } 777 }
769 778
770 G4Exception( fMethod, "GeomNav0003", excepti 779 G4Exception( fMethod, "GeomNav0003", exceptionType, msg);
771 } 780 }
772 781
773 namespace G4NavigationLogger_Namespace << 782 namespace G4NavigationLogger_Namespace {
774 { <<
775 const G4String EInsideNames[3] = { "kOutside 783 const G4String EInsideNames[3] = { "kOutside", "kSurface", "kInside" };
776 } 784 }
777 785
778 void G4NavigationLogger:: << 786 void
779 ReportVolumeAndIntersection( std::ostream& os, << 787 G4NavigationLogger::ReportVolumeAndIntersection( std::ostream& os,
780 const G4ThreeVect << 788 const G4ThreeVector& localPoint,
781 const G4ThreeVect << 789 const G4ThreeVector& localDirection,
782 const G4VPhysical << 790 const G4VPhysicalVolume* physical ) const
783 { << 791 {
784 G4String fMethod = fId + "::ComputeStep()"; << 792 G4String fMethod = fId + "::ComputeStep() / output from G4NavigatorLogger::ReportVolumeAndIntersection";
785 const G4LogicalVolume* logicalVol = physical << 793 const G4LogicalVolume* logicalVol = physical
786 ? physical << 794 ? physical->GetLogicalVolume() : 0;
787 const G4VSolid* solid = logicalVol != nullpt << 795 const G4VSolid* solid = logicalVol
788 ? logicalVol->GetSolid << 796 ? logicalVol->GetSolid() : 0;
789 if( solid == nullptr ) << 797 if( solid == 0 )
790 { 798 {
791 os << " ERROR> Solid is not available. Lo << 799 os << " ERROR> Solid is not available. Logical Volume = " << logicalVol << std::endl;
792 << logicalVol << std::endl; <<
793 return; 800 return;
794 } 801 }
795 const G4double kCarTolerance = solid->GetTol 802 const G4double kCarTolerance = solid->GetTolerance();
796 803
797 // Double check reply - it should be kInfini 804 // Double check reply - it should be kInfinity
798 const G4double distToOut = solid->DistanceTo 805 const G4double distToOut = solid->DistanceToOut(localPoint, localDirection);
799 const G4double distToOutNeg = solid->Distanc << 806 const G4double distToOutNeg = solid->DistanceToOut(localPoint, -localDirection);
800 << 807
801 const EInside inSolid = solid->Inside(loc 808 const EInside inSolid = solid->Inside(localPoint);
802 const G4double safetyToIn = solid->Distance 809 const G4double safetyToIn = solid->DistanceToIn(localPoint);
803 const G4double safetyToOut = solid->Distance 810 const G4double safetyToOut = solid->DistanceToOut(localPoint);
804 811
805 const G4double distToInPos = solid->Distance 812 const G4double distToInPos = solid->DistanceToIn(localPoint, localDirection);
806 const G4double distToInNeg = solid->Distance 813 const G4double distToInNeg = solid->DistanceToIn(localPoint, -localDirection);
807 814
808 const G4ThreeVector exitNormal = solid->Surf 815 const G4ThreeVector exitNormal = solid->SurfaceNormal(localPoint);
809 816
810 // Double check whether points nearby are in 817 // Double check whether points nearby are in/surface/out
811 const G4double epsilonDist = 1000.0 * kCarTo 818 const G4double epsilonDist = 1000.0 * kCarTolerance;
812 const G4ThreeVector PointPlusDir = localPoin << 819 const G4ThreeVector PointPlusDir = localPoint + epsilonDist * localDirection;
813 const G4ThreeVector PointMinusDir = localPoi 820 const G4ThreeVector PointMinusDir = localPoint - epsilonDist * localDirection;
814 const G4ThreeVector PointPlusNorm = localPoi << 821 const G4ThreeVector PointPlusNorm = localPoint + epsilonDist * exitNormal;
815 const G4ThreeVector PointMinusNorm = localPo 822 const G4ThreeVector PointMinusNorm = localPoint - epsilonDist * exitNormal;
816 823
817 const EInside inPlusDir = solid->Inside(Poin << 824 const EInside inPlusDir= solid->Inside(PointPlusDir);
818 const EInside inMinusDir = solid->Inside(Poi << 825 const EInside inMinusDir= solid->Inside(PointMinusDir);
819 const EInside inPlusNorm = solid->Inside(Poi << 826 const EInside inPlusNorm= solid->Inside(PointPlusNorm);
820 const EInside inMinusNorm = solid->Inside(Po << 827 const EInside inMinusNorm= solid->Inside(PointMinusNorm);
821 828
822 // Basic information 829 // Basic information
823 os << " Current physical volume = " << phy 830 os << " Current physical volume = " << physical->GetName() << G4endl;
824 os << " Position (loc) = " << localPoint 831 os << " Position (loc) = " << localPoint << G4endl
825 << " Direction (dir) = " << localDirect << 832 << " Direction (dir) = " << localDirection << G4endl;
826 os << " For confirmation:" << G4endl; 833 os << " For confirmation:" << G4endl;
827 os << " Response of DistanceToOut (loc, +d 834 os << " Response of DistanceToOut (loc, +dir)= " << distToOut << G4endl;
828 os << " Response of DistanceToOut (loc, -d 835 os << " Response of DistanceToOut (loc, -dir)= " << distToOutNeg << G4endl;
829 836
830 os << " Inside responds = " << inSolid << 837 os << " Inside responds = " << inSolid << " , ie: ";
831 if( inSolid == kOutside ) << 838 if( inSolid == kOutside ) {
832 { <<
833 os << " Outside -- a problem, as observed 839 os << " Outside -- a problem, as observed in " << fMethod << G4endl;
834 } << 840 } else if( inSolid == kSurface ) {
835 else if( inSolid == kSurface ) <<
836 { <<
837 os << " Surface -- unexpected / inconsiste 841 os << " Surface -- unexpected / inconsistent response ! " << G4endl;
838 } << 842 } else {
839 else <<
840 { <<
841 os << " Inside -- unexpected / inconsiste 843 os << " Inside -- unexpected / inconsistent response ! " << G4endl;
842 } 844 }
843 os << " Obtain safety(ToIn) = " << safetyT 845 os << " Obtain safety(ToIn) = " << safetyToIn << G4endl;
844 os << " Obtain safety(ToOut) = " << safety 846 os << " Obtain safety(ToOut) = " << safetyToOut << G4endl;
845 os << " Response of DistanceToIn (loc, +dir) 847 os << " Response of DistanceToIn (loc, +dir)= " << distToInPos << G4endl;
846 os << " Response of DistanceToIn (loc, -dir) 848 os << " Response of DistanceToIn (loc, -dir)= " << distToInNeg << G4endl;
847 849
848 os << " Exit Normal at loc = " << exitNormal 850 os << " Exit Normal at loc = " << exitNormal << G4endl;
849 os << " Dir . Normal = " << exitNormal 851 os << " Dir . Normal = " << exitNormal.dot( localDirection );
850 os << G4endl; 852 os << G4endl;
851 853
852 os << " Checking points moved from position << 854 os << " Checking points moved from position by distance/dir. Solid responses: " << G4endl
853 << " Solid responses: " << G4endl << 855 << " +eps in direction : " << G4NavigationLogger_Namespace::EInsideNames[inPlusDir]
854 << " +eps in direction : " << 856 << " +eps in Normal : " << G4NavigationLogger_Namespace::EInsideNames[inPlusNorm] << G4endl
855 << G4NavigationLogger_Namespace::EInsideN << 857 << " -eps in direction : " << G4NavigationLogger_Namespace::EInsideNames[inMinusDir]
856 << " +eps in Normal : " << 858 << " -eps in Normal : " << G4NavigationLogger_Namespace::EInsideNames[inMinusNorm] << G4endl;
857 << G4NavigationLogger_Namespace::EInsideN <<
858 << " -eps in direction : " <<
859 << G4NavigationLogger_Namespace::EInsideN <<
860 << " -eps in Normal : " <<
861 << G4NavigationLogger_Namespace::EInsideN <<
862 859
863 os << " Parameters of solid: " << G4endl 860 os << " Parameters of solid: " << G4endl;
864 os << *solid; 861 os << *solid;
865 os << "===================================== 862 os << "============================================================";
866 } 863 }
867 864