Changeset 991 for trunk/source/processes/electromagnetic/standard/src/G4WentzelVIModel.cc

Timestamp:

Apr 17, 2009, 12:17:14 PM (15 years ago)

Author:

garnier

Message:

update

File:

• trunk/source/processes/electromagnetic/standard/src/G4WentzelVIModel.cc (modified) (5 diffs)

trunk/source/processes/electromagnetic/standard/src/G4WentzelVIModel.cc

@@ -24 +24 @@
 // ********************************************************************
 //
-// $Id: G4WentzelVIModel.cc,v 1.17 2009/02/19 19:17:15 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+// $Id: G4WentzelVIModel.cc,v 1.16 2008/11/19 11:47:50 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-02 $
 //
 // -------------------------------------------------------------------
@@ -60 +60 @@
 #include "G4LossTableManager.hh"
 #include "G4ParticleChangeForMSC.hh"
-//#include "G4TransportationManager.hh"
-//#include "G4SafetyHelper.hh"
+#include "G4TransportationManager.hh"
+#include "G4SafetyHelper.hh"
 #include "G4PhysicsTableHelper.hh"
 #include "G4ElementVector.hh"
@@ -138 +138 @@
       fParticleChange = new G4ParticleChangeForMSC();
 
-    InitialiseSafetyHelper();
+    safetyHelper = G4TransportationManager::GetTransportationManager()
+      ->GetSafetyHelper();
+    safetyHelper->InitialiseHelper();
   }
 }
@@ -275 +277 @@
   // i.e. when it is needed for optimization purposes
   if(stepStatus != fGeomBoundary && presafety < tlimitminfix) 
-    presafety = ComputeSafety(sp->GetPosition(), tlimit); 
+    presafety = safetyHelper->ComputeSafety(sp->GetPosition()); 
   /*
   G4cout << "G4WentzelVIModel::ComputeTruePathLengthLimit tlimit= " 
@@ -568 +570 @@
 
     if(r > tlimitminfix) {
-      pos /= r;
-      ComputeDisplacement(fParticleChange, pos, r, safety);
+      G4ThreeVector Position = *(fParticleChange->GetProposedPosition());
+      G4double fac= 1.;
+      if(r >= safety) {
+        //  ******* so safety is computed at boundary too ************
+        G4double newsafety = 
+          safetyHelper->ComputeSafety(Position) - tlimitminfix;
+        if(newsafety <= 0.0) fac = 0.0;
+        else if(r > newsafety) fac = newsafety/r ;
+        //G4cout << "NewSafety= " << newsafety << " fac= " << fac 
+        // << " r= " << r << " sint= " << sint << " pos " << Position << G4endl;
+      }  
+
+      if(fac > 0.) {
+        // compute new endpoint of the Step
+        G4ThreeVector newPosition = Position + fac*pos;
+
+        // check safety after displacement
+        G4double postsafety = safetyHelper->ComputeSafety(newPosition);
+
+        // displacement to boundary
+        if(postsafety <= 0.0) {
+          safetyHelper->Locate(newPosition, newDirection);
+
+          // not on the boundary
+        } else { 
+          safetyHelper->ReLocateWithinVolume(newPosition);
+          // if(fac < 1.0) G4cout << "NewPosition " << newPosition << G4endl;
+        }
+     
+        fParticleChange->ProposePosition(newPosition);
+      } 
     }
   }