4 ////////////////////////////////////////////////////////////////////
6 // WrapNorml.hh - Sara Vanini
8 // Wrapper for computing normal unit-vector in global coordinates.
10 // Fluka requires normal vector exiting from final position (of
11 // particle) volume, that is: entering in volume of initial position.
12 // Geant4 always computes the normal vector exiting from the volume.
13 // In GetLocalExitNormal() call the volume is the pre-step volume
14 // (so G4 normal vector sign is opposite of fluka-required normal).
15 // If IsLocalExitNormalValid=false, normal value is computed from
16 // init-step volume (in this case sign must change), or from
17 // end-step volume (the sign is the same). Normal vector is computed
18 // always on boundary between volumes, in global coordinates (to take
19 // rotation of parameterised volumes in hierarchy in consideration).
20 // So: nrmlwr returns inwards pointing unit normal of the shape for
21 // surface closest to the point returned by the navigator (last step
26 // modified 7/VI/00 for boundary-crossing in case of relocation
27 // modified 5/VII/00 geometry error on boundary fixed
28 // modified 24.10.01: by I. Hrivnacova
29 // functions declarations separated from implementation
30 // (moved to Wrappers.hh);
32 ////////////////////////////////////////////////////////////////////
34 #include "Wrappers.hh"
35 #include "WrapUtils.hh"
36 #include "FGeometryInit.hh"
37 #include "G4VPhysicalVolume.hh"
38 #include "FluggNavigator.hh"
39 #include "G4ThreeVector.hh"
43 void nrmlwr(G4double& /*pSx*/, G4double& /*pSy*/, G4double& /*pSz*/,
44 G4double& pVx, G4double& pVy, G4double& pVz,
45 G4double* norml, const G4int& oldReg,
46 const G4int& /*newReg*/, G4int& flagErr)
49 #ifdef G4GEOMETRY_DEBUG
50 G4cout << "============ NRMLWR-DBG =============" << G4endl;
57 static FGeometryInit * ptrGeoInit;
58 ptrGeoInit = FGeometryInit::GetInstance();
59 FluggNavigator * ptrNavig = ptrGeoInit->getNavigatorForTracking();
62 G4ThreeVector normalLoc;
63 G4ThreeVector normalGlob;
67 normalLoc=ptrNavig->GetLocalExitNormal(&valid);
71 //global cooordinates normal
72 normalGlob = ptrNavig->GetLocalToGlobalTransform().
73 TransformAxis(normalLoc);
76 G4VPhysicalVolume *touchvolume;
77 G4ThreeVector theLocalPoint;
79 if(ptrNavig->EnteredDaughterVolume()) {
80 //volume from touchable history
81 G4TouchableHistory *ptrTouchableHistory=ptrGeoInit->
82 GetTouchableHistory();
83 touchvolume = ptrTouchableHistory->GetVolume();
85 //local point from navigator and normal
86 theLocalPoint = ptrNavig->GetCurrentLocalCoordinate();
87 normalLoc = touchvolume->GetLogicalVolume()->GetSolid()->
88 SurfaceNormal(theLocalPoint);
90 //global cooordinates normal
91 normalGlob = ptrNavig->GetLocalToGlobalTransform().
92 TransformAxis(normalLoc);
95 //volume from old history
96 const G4NavigationHistory * ptrOldNavHist =
97 ptrGeoInit->GetOldNavHist()->GetHistory();
98 touchvolume = ptrOldNavHist->GetTopVolume();
101 // Old history has been reseted by LOOKZ relocation,
102 // so is necessary to track back and forward to find
103 // the right histories.
106 //////////// COMPUTE STEP BACKWARD //////////////////
107 G4ThreeVector theGlobalPoint = ptrNavig->GetLocalToGlobalTransform().
108 TransformPoint(ptrNavig->GetCurrentLocalCoordinate());
110 //compute step and new location
111 G4ThreeVector pVec(pVx,pVy,pVz);
113 G4double appStep = 0;
115 G4bool onBoundary = false;
116 G4double physStep = 1000000000;
118 G4ThreeVector partLoc = theGlobalPoint;
121 #ifdef G4GEOMETRY_DEBUG
122 G4cout << "Old history not found" << G4endl;
123 G4cout << "* NRML needs boundary-crossing: computing step backward..."
127 //compute step and location
128 newRegStep=StepAndLocation(partLoc,pVec,physStep,
129 appStep,safety,onBoundary,
132 G4PhysicalVolumeStore * pVolStore = G4PhysicalVolumeStore::
135 if(appStep<physStep && newRegStep!=G4int(pVolStore->size())+1) {
136 //end-step point is on boundary between volumes;
137 //==> update step-histories
139 #ifdef G4GEOMETRY_DEBUG
140 G4cout << "* updating step-histories" << G4endl;
143 G4TouchableHistory * ptrTouchHist =
144 ptrGeoInit->GetTouchableHistory();
145 ptrGeoInit->UpdateHistories(ptrTouchHist->GetHistory(),2);
148 #ifdef G4GEOMETRY_DEBUG
149 G4cout << "ERROR! Boundary not found" << G4endl;
153 #ifdef G4GEOMETRY_DEBUG
154 G4cout << "* computing step forward..." << G4endl;
160 //compute step and location for boundary crossing
161 newRegStep=StepAndLocation(partLoc,pVec,physStep,
162 appStep,safety,onBoundary,fErr,oldReg);
163 if(appStep<physStep) {
164 //end-step point is on boundary between volumes;
165 //==> update step-histories
167 #ifdef G4GEOMETRY_DEBUG
168 G4cout << "* updating step-histories" << G4endl;
171 G4TouchableHistory * ptrTouchHist =
172 ptrGeoInit->GetTouchableHistory();
173 ptrGeoInit->UpdateHistories(ptrTouchHist->GetHistory(),2);
177 // now touchvolume exist.
178 // local point from navigator and global point
179 // N.B. if particle has exited world volume,
180 // FluggNavigator doesn't update lastLocatedPoint.
181 // So be carefull in building geometry always to have a
182 // big world volume that fluka won't exit.
184 touchvolume = ptrOldNavHist->GetTopVolume();
186 G4ThreeVector theGlobalPoint = ptrNavig->GetLocalToGlobalTransform().
187 TransformPoint(ptrNavig->GetCurrentLocalCoordinate());
188 theLocalPoint = ptrOldNavHist->GetTopTransform().
189 TransformPoint(theGlobalPoint);
190 normalLoc = (touchvolume->GetLogicalVolume()->GetSolid()->
191 SurfaceNormal(theLocalPoint));
194 //global cooordinates normal
195 normalGlob = ptrOldNavHist->GetTopTransform().
196 Inverse().TransformAxis(normalLoc);
201 norml[0]=G4double(normalGlob.x());
202 norml[1]=G4double(normalGlob.y());
203 norml[2]=G4double(normalGlob.z());
206 #ifdef G4GEOMETRY_DEBUG
207 G4cout << "Normal: " << normalGlob << G4endl;