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