Index: trunk/source/processes/hadronic/management/src/G4HadronicProcess.cc
===================================================================
--- trunk/source/processes/hadronic/management/src/G4HadronicProcess.cc	(revision 1228)
+++ trunk/source/processes/hadronic/management/src/G4HadronicProcess.cc	(revision 1315)
@@ -45,16 +45,15 @@
 #include "G4StableIsotopes.hh"
 #include "G4HadTmpUtil.hh"
-
-#include "G4HadLeadBias.hh"
+#include "G4NucleiProperties.hh"
+
+//VI #include "G4HadLeadBias.hh"
 #include "G4HadronicException.hh"
-#include "G4HadReentrentException.hh"
-
-#include "G4HadronicWhiteBoard.hh"
-
-#include "G4HadSignalHandler.hh"
-
+//VI #include "G4HadReentrentException.hh"
+//VI #include "G4HadronicWhiteBoard.hh"
+//VI #include "G4HadSignalHandler.hh"
 #include "G4HadronicProcessStore.hh"
 
-#include <typeinfo>
+
+/*VI
 
 namespace G4HadronicProcess_local
@@ -65,4 +64,5 @@
   }
 } 
+VI*/
 
 G4IsoParticleChange * G4HadronicProcess::theIsoResult = 0;
@@ -78,11 +78,12 @@
 //////////////////////////////////////////////////////////////////
 
-G4HadronicProcess::G4HadronicProcess( const G4String &processName,
-                                      G4ProcessType   aType ) :
-G4VDiscreteProcess( processName, aType)
-{ 
+G4HadronicProcess::G4HadronicProcess(const G4String& processName,
+                                     G4ProcessType aType)
+ :G4VDiscreteProcess(processName, aType)
+{
   ModelingState = 0;
   isoIsOnAnyway = -1;
   theTotalResult = new G4ParticleChange();
+  theInteraction = 0;
   theCrossSectionDataStore = new G4CrossSectionDataStore();
   G4HadronicProcessStore::Instance()->Register(this);
@@ -90,5 +91,10 @@
   xBiasOn = false;
   G4HadronicProcess_debug_flag = false;
-  if(getenv("SwitchLeadBiasOn")) theBias.push_back(new G4HadLeadBias());
+  epReportLevel = 0;
+  epCheckLevels.first = DBL_MAX;
+  epCheckLevels.second = DBL_MAX;
+  levelsSetByProcess = false;
+
+  //VI   if(getenv("SwitchLeadBiasOn")) theBias.push_back(new G4HadLeadBias());
 }
 
@@ -100,5 +106,5 @@
   std::for_each(theProductionModels.begin(),
                 theProductionModels.end(), G4Delete());
-  std::for_each(theBias.begin(), theBias.end(), G4Delete());
+  //VI  std::for_each(theBias.begin(), theBias.end(), G4Delete());
  
   delete theOldIsoResult; 
@@ -109,8 +115,10 @@
 void G4HadronicProcess::RegisterMe( G4HadronicInteraction *a )
 { 
+  if(!a) { return; }
   try{GetManagerPointer()->RegisterMe( a );}   
   catch(G4HadronicException & aE)
   {
-    aE.Report(std::cout);
+    G4cout << "Unrecoverable error in " << GetProcessName() 
+	   << " to register " << a->GetModelName() << G4endl;
     G4Exception("G4HadronicProcess", "007", FatalException,
     "Could not register G4HadronicInteraction");
@@ -136,5 +144,5 @@
   try
   {
-    ModelingState = 1;    
+    //VI ModelingState = 1;    
     theLastCrossSection = aScaleFactor* 
       theCrossSectionDataStore->GetCrossSection(aTrack.GetDynamicParticle(), 
@@ -143,10 +151,10 @@
   catch(G4HadronicException aR)
   { 
-    aR.Report(G4cout); 
+    DumpState(aTrack,"GetMeanFreePath");
     G4Exception("G4HadronicProcess", "007", FatalException,
     "G4HadronicProcess::GetMeanFreePath failed");
   } 
   G4double res = DBL_MAX;
-  if( theLastCrossSection > 0.0 ) res = 1.0/theLastCrossSection;
+  if( theLastCrossSection > 0.0 ) { res = 1.0/theLastCrossSection; }
   return res;
 }
@@ -165,10 +173,10 @@
 {
   // Debugging stuff
-  if(G4HadronicProcess_debug_flag) 
-         std::cout << "@@@@ hadronic process start "<< std::endl;
+  //VI if(G4HadronicProcess_debug_flag) 
+  //VI       std::cout << "@@@@ hadronic process start "<< std::endl;
   // G4cout << theNumberOfInteractionLengthLeft<<G4endl;
-#ifndef G4HadSignalHandler_off
-  G4HadSignalHandler aHandler(G4HadronicProcess_local::G4HadronicProcessHandler_1);
-#endif
+  //VI #ifndef G4HadSignalHandler_off
+  //VI G4HadSignalHandler aHandler(G4HadronicProcess_local::G4HadronicProcessHandler_1);
+  //VI #endif
 
   if(aTrack.GetTrackStatus() != fAlive && aTrack.GetTrackStatus() != fSuspend) {
@@ -176,7 +184,8 @@
         aTrack.GetTrackStatus() == fKillTrackAndSecondaries ||
         aTrack.GetTrackStatus() == fPostponeToNextEvent) {
-      G4cerr << "G4HadronicProcess: track in unusable state - "
+      G4cout << "G4HadronicProcess: track in unusable state - "
              << aTrack.GetTrackStatus() << G4endl;
-      G4cerr << "G4HadronicProcess: returning unchanged track " << G4endl;
+      G4cout << "G4HadronicProcess: returning unchanged track " << G4endl;
+      DumpState(aTrack,"PostStepDoIt");
       G4Exception("G4HadronicProcess", "001", JustWarning, "bailing out");
     }
@@ -224,10 +233,10 @@
   catch(G4HadronicException & aR)
   {
-    aR.Report(G4cout);
-    G4cout << "Unrecoverable error for:"<<G4endl;
-    G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
-    G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
-    G4cout << " - Particle type = "
-    <<aParticle->GetDefinition()->GetParticleName()<<G4endl;
+    DumpState(aTrack,"SampleZandA");
+    //VI    G4cout << "Unrecoverable error for:"<<G4endl;
+    //VI G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
+    //VI G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
+    //VI G4cout << " - Particle type = "
+    //VI <<aParticle->GetDefinition()->GetParticleName()<<G4endl;
     G4Exception("G4HadronicProcess", "007", FatalException,
     "PostStepDoIt failed on element selection.");
@@ -236,15 +245,16 @@
   try
   {
-    theInteraction = ChooseHadronicInteraction( kineticEnergy,
-    aMaterial, anElement );
+    theInteraction = 
+      ChooseHadronicInteraction( kineticEnergy, aMaterial, anElement );
   }
   catch(G4HadronicException & aE)
   {
-    aE.Report(std::cout);
-    G4cout << "Unrecoverable error for:"<<G4endl;
-    G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
-    G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
-    G4cout << " - Particle type = " 
-           << aParticle->GetDefinition()->GetParticleName()<<G4endl;
+    G4cout << "Target element "<<anElement->GetName()<<"  Z= " 
+	   << targetNucleus.GetZ() << "  A= " << targetNucleus.GetN() << G4endl;
+    DumpState(aTrack,"ChooseHadronicInteraction");
+    //VI G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
+    //VI G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
+    //VI G4cout << " - Particle type = " 
+    //VI       << aParticle->GetDefinition()->GetParticleName()<<G4endl;
     G4Exception("G4HadronicProcess", "007", FatalException,
     "ChooseHadronicInteraction failed.");
@@ -264,5 +274,7 @@
       // Call the interaction
       result = theInteraction->ApplyYourself( thePro, targetNucleus);
-    }
+      ++reentryCount;
+    }
+    /*VI
     catch(G4HadReentrentException aR)
     {
@@ -281,50 +293,56 @@
              << aParticle->GetDefinition()->GetParticleName() << G4endl;
       result = 0; // here would still be leaking...
-      if(reentryCount>100)
-      {
-        G4Exception("G4HadronicProcess", "007", FatalException,
-        "GetHadronicProcess: Reentering ApplyYourself too often - PostStepDoIt failed.");  
-      }
       G4Exception("G4HadronicProcess", "007", FatalException,
       "GetHadronicProcess: PostStepDoIt failed (Reentering ApplyYourself not yet supported.)");  
     }
+    VI*/
     catch(G4HadronicException aR)
     {
-      G4HadronicWhiteBoard & theBoard = G4HadronicWhiteBoard::Instance();
-      theBoard.SetProjectile(thePro);
-      theBoard.SetTargetNucleus(targetNucleus);
-      theBoard.SetProcessName(GetProcessName());
-      theBoard.SetModelName(theInteraction->GetModelName());
-
-      aR.Report(G4cout);
-      G4cout << " G4HadronicProcess failed in ApplyYourself call for" 
-             << G4endl;
-      G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
-      G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
-      G4cout << " - Particle type = " 
-             << aParticle->GetDefinition()->GetParticleName() << G4endl;
+      //VI G4HadronicWhiteBoard & theBoard = G4HadronicWhiteBoard::Instance();
+      //VI theBoard.SetProjectile(thePro);
+      //VI theBoard.SetTargetNucleus(targetNucleus);
+      //VI theBoard.SetProcessName(GetProcessName());
+      //VI theBoard.SetModelName(theInteraction->GetModelName());
+
+      G4cout << "Call for " << theInteraction->GetModelName() << G4endl;
+      G4cout << "Target element "<<anElement->GetName()<<"  Z= " 
+	     << targetNucleus.GetZ() << "  A= " << targetNucleus.GetN() << G4endl;
+      DumpState(aTrack,"ApplyYourself");
+      //VI G4cout << " - Particle energy[GeV] = "<< originalEnergy/GeV<<G4endl;
+      //VI G4cout << " - Material = "<<aMaterial->GetName()<<G4endl;
+      //VI G4cout << " - Particle type = " 
+      //VI        << aParticle->GetDefinition()->GetParticleName() << G4endl;
       G4Exception("G4HadronicProcess", "007", FatalException,
       "PostStepDoIt failed.");
     }
+    if(reentryCount>100) {
+      G4cout << "Call for " << theInteraction->GetModelName() << G4endl;
+      G4cout << "Target element "<<anElement->GetName()<<"  Z= " 
+	     << targetNucleus.GetZ() << "  A= " << targetNucleus.GetN() << G4endl;
+      DumpState(aTrack,"ApplyYourself");
+      G4Exception("G4HadronicProcess", "007", FatalException,
+		  "Reentering ApplyYourself too often - PostStepDoIt failed.");  
+    }
   }
   while(!result);
 
-  if(!ModelingState && !getenv("BypassAllSafetyChecks") ) 
-  {
-    G4cout << "ERROR IN EXECUTION -- HADRONIC PROCESS STATE NOT VALID"<<G4endl;
-    G4cout << "Result will be of undefined quality."<<G4endl;
-  }
+  //VI if(!ModelingState && !getenv("BypassAllSafetyChecks") ) 
+  //VI{
+  //VI  G4cout << "ERROR IN EXECUTION -- HADRONIC PROCESS STATE NOT VALID"<<G4endl;
+  //VI  G4cout << "Result will be of undefined quality."<<G4endl;
+  //VI}
 
   // NOT USED ?? Projectile particle has changed character during interaction
-  if(result->GetStatusChange() == isAlive && 
-     thePro.GetDefinition() != aTrack.GetDefinition())
-  {
-    G4DynamicParticle * aP = 
-             const_cast<G4DynamicParticle *>(aTrack.GetDynamicParticle());
-    aP->SetDefinition(const_cast<G4ParticleDefinition *>(thePro.GetDefinition()));
-  }
+  //VI if(result->GetStatusChange() == isAlive && 
+  //VI    thePro.GetDefinition() != aTrack.GetDefinition())
+  //VI {
+  //VI  G4DynamicParticle * aP = 
+  //VI           const_cast<G4DynamicParticle *>(aTrack.GetDynamicParticle());
+  //VI   aP->SetDefinition(const_cast<G4ParticleDefinition *>(thePro.GetDefinition()));
+  //VI }
 
   result->SetTrafoToLab(thePro.GetTrafoToLab());
 
+  /*VI
   if(getenv("HadronicDoitLogging") )
   {
@@ -338,4 +356,5 @@
 	   << G4endl;
   }
+  VI*/
 
   ClearNumberOfInteractionLengthLeft();
@@ -348,4 +367,5 @@
   }
   
+  /*VI
   G4double e=aTrack.GetKineticEnergy();
   ModelingState = 0;
@@ -357,14 +377,17 @@
     }
   }
-
+  VI*/
   // Put hadronic final state particles into G4ParticleChange
 
   FillTotalResult(result, aTrack);
-  if(G4HadronicProcess_debug_flag) 
-    std::cout << "@@@@ hadronic process end "<< std::endl;
-    
+  //VI if(G4HadronicProcess_debug_flag) 
+  //VI  std::cout << "@@@@ hadronic process end "<< std::endl;
+
+  if (epReportLevel > 0) CheckEnergyMomentumConservation(aTrack, targetNucleus);
+ 
   return theTotalResult;
 }
 
+#include <typeinfo>
 
 G4HadFinalState* 
@@ -487,5 +510,5 @@
   theTotalResult->SetSecondaryWeightByProcess(true);
   theTotalResult->ProposeTrackStatus(fAlive);
-  G4double rotation = 2.*pi*G4UniformRand();
+  G4double rotation = CLHEP::twopi*G4UniformRand();
   G4ThreeVector it(0., 0., 1.);
 
@@ -520,4 +543,5 @@
       if(xBiasOn)
       {
+        DumpState(aT,"FillTotalResult");
         G4Exception("G4HadronicProcess", "007", FatalException,
         "Cannot cross-section bias a process that suspends tracks.");
@@ -575,5 +599,8 @@
   else
   {
-    G4cerr << "Track status is "<< aR->GetStatusChange()<<G4endl;
+    G4cout << "Call for " << theInteraction->GetModelName() << G4endl;
+    G4cout << "Target Z= " 
+	   << targetNucleus.GetZ() << "  A= " << targetNucleus.GetN() << G4endl;
+    DumpState(aT,"FillTotalResult");
     G4Exception("G4HadronicProcess", "007", FatalException,
     "use of unsupported track-status.");
@@ -691,13 +718,108 @@
       (it != "PositronNuclear") )
     {
-      G4Exception("G4HadronicProcess", "007", FatalException,
+      G4Exception("G4HadronicProcess::BiasCrossSectionByFactor", "007", FatalException,
 		  "Cross-section biasing available only for gamma and electro nuclear reactions.");
     }
   if(aScale<100)
     {
-      G4Exception("G4HadronicProcess", "001", JustWarning,
+      G4Exception("G4HadronicProcess::BiasCrossSectionByFactor", "001", JustWarning,
 		  "Cross-section bias readjusted to be above safe limit. New value is 100");
       aScaleFactor = 100.;
     }
 }
+
+void 
+G4HadronicProcess::CheckEnergyMomentumConservation(const G4Track& aTrack,
+                                                   const G4Nucleus& aNucleus)
+{
+  G4double targetMass = G4NucleiProperties::GetNuclearMass(aNucleus.GetN(), aNucleus.GetZ());
+  G4LorentzVector projectile4mom = aTrack.GetDynamicParticle()->Get4Momentum();
+  G4LorentzVector target4mom(0, 0, 0, targetMass);
+  G4LorentzVector initial4mom = projectile4mom + target4mom;
+
+  G4Track* sec;
+  G4LorentzVector final4mom;
+  G4int nSec = theTotalResult->GetNumberOfSecondaries();
+  for (G4int i = 0; i < nSec; i++) {
+    sec = theTotalResult->GetSecondary(i);
+    final4mom += sec->GetDynamicParticle()->Get4Momentum();
+  }
+
+  G4LorentzVector diff = initial4mom - final4mom;
+  G4double absolute = diff.e();
+  G4double relative = absolute/aTrack.GetKineticEnergy();
+
+  G4String processName = GetProcessName();
+  G4HadronicInteraction* theModel = GetHadronicInteraction();
+  G4String modelName("none");
+  if (theModel) modelName = theModel->GetModelName();
+
+  std::pair<G4double, G4double> checkLevels = epCheckLevels;;
+  if (!levelsSetByProcess) {
+    if (theModel) checkLevels = theModel->GetEnergyMomentumCheckLevels();
+  }
+
+  G4bool relPass = false;
+  G4String relResult = "fail";
+  if (std::abs(relative) < checkLevels.first) {
+    relPass = true;
+    relResult = "pass";
+  }
+
+  G4bool absPass = false;
+  G4String absResult = "fail";
+  if (std::abs(absolute) < checkLevels.second) {
+    absPass = true;
+    absResult = "pass";
+  }
+
+  // Options for level of reporting detail:
+  //  0. off
+  //  1. report only when E/p not conserved
+  //  2. report regardless of E/p conservation
+  //  3. report only when E/p not conserved, with model names, process names, and limits 
+  //  4. report regardless of E/p conservation, with model names, process names, and limits
+
+  if(epReportLevel == 4) {
+    G4cout << " Process: " << processName << " , Model: " <<  modelName << G4endl; 
+    G4cout << " relative limit " << checkLevels.first << " relative value = "
+           << relative << " " << relResult << G4endl;
+    G4cout << " absolute limit " << checkLevels.second << " absolute value = "
+           << absolute << " " << absResult << G4endl;
+
+  } else if(epReportLevel == 3) {
+    if (!absPass || !relPass) {
+      G4cout << " Process: " << processName << " , Model: " <<  modelName << G4endl; 
+      G4cout << " relative limit " << checkLevels.first << " relative value = "
+             << relative << " " << relResult << G4endl;  
+      G4cout << " absolute limit " << checkLevels.second << " absolute value = "
+             << absolute << " " << absResult << G4endl;
+    }
+
+  } else if(epReportLevel == 2) {
+    G4cout << " relative value = " << relative << " " << relPass
+           << " absolute value = " << absolute << " " << absPass << G4endl;
+
+  } else if(epReportLevel == 1) {
+    if (!absPass || !relPass) {
+      G4cout << " relative value = " << relative << " " << relPass
+             << " absolute value = " << absolute << " " << absPass << G4endl;
+    }
+  }
+}
+
+
+void G4HadronicProcess::DumpState(const G4Track& aTrack, const G4String& method)
+{
+  G4cout << "Unrecoverable error in the method " << method << " of " 
+	 << GetProcessName() << G4endl;
+  G4cout << "TrackID= "<< aTrack.GetTrackID() << "  ParentID= " << aTrack.GetParentID()
+	 << "  " << aTrack.GetDefinition()->GetParticleName() << G4endl;
+  G4cout << "Ekin(GeV)= " << aTrack.GetKineticEnergy()/CLHEP::GeV 
+	 << ";  direction= " << aTrack.GetMomentumDirection() << G4endl; 
+  G4cout << "Position(mm)= " << aTrack.GetPosition()/CLHEP::mm 
+	 << ";  material " << aTrack.GetMaterial()->GetName() << G4endl; 
+  G4cout << "PhysicalVolume  <" << aTrack.GetVolume()->GetName() << ">" << G4endl;
+} 
+
 /* end of file */
Index: trunk/source/processes/hadronic/management/src/G4HadronicProcessStore.cc
===================================================================
--- trunk/source/processes/hadronic/management/src/G4HadronicProcessStore.cc	(revision 1228)
+++ trunk/source/processes/hadronic/management/src/G4HadronicProcessStore.cc	(revision 1315)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4HadronicProcessStore.cc,v 1.13 2009/12/02 15:56:02 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-03 $
+// $Id: G4HadronicProcessStore.cc,v 1.14 2010/04/21 17:55:25 dennis Exp $
+// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
 //
 // -------------------------------------------------------------------
@@ -57,4 +57,5 @@
 #include "G4HadronicInteractionRegistry.hh"
 #include "G4CrossSectionDataSetRegistry.hh"
+#include "G4HadronicEPTestMessenger.hh"
 
 G4HadronicProcessStore* G4HadronicProcessStore::theInstance = 0;
@@ -78,4 +79,5 @@
   G4HadronicInteractionRegistry::Instance()->Clean();
   G4CrossSectionDataSetRegistry::Instance()->Clean();
+  delete theEPTestMessenger;
 }
 
@@ -127,4 +129,5 @@
   verbose = 1;
   buildTableStart = true;
+  theEPTestMessenger = new G4HadronicEPTestMessenger(this);
 }
 
@@ -367,5 +370,5 @@
     }
   }
-  //G4cout << "G4HadronicProcessStore::Register hadronic " << n_proc
+  //  G4cout << "G4HadronicProcessStore::Register hadronic " << n_proc
   //	 << "  " << proc->GetProcessName() << G4endl;
   n_proc++;
@@ -664,3 +667,37 @@
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
+
+void G4HadronicProcessStore::SetEpReportLevel(G4int level)
+{
+  G4cout << " Setting energy/momentum report level to " << level 
+         << " for " << process.size() << " hadronic processes " << G4endl;
+  for (G4int i = 0; i < G4int(process.size()); i++) {
+    process[i]->SetEpReportLevel(level);
+  }
+}
+
+
+void G4HadronicProcessStore::SetProcessAbsLevel(G4double abslevel)
+{
+  G4cout << " Setting absolute energy/momentum test level to " << abslevel << G4endl;
+  G4double rellevel = 0.0;
+  G4HadronicProcess* theProcess = 0;
+  for (G4int i = 0; i < G4int(process.size()); i++) {
+    theProcess = process[i];
+    rellevel = theProcess->GetEnergyMomentumCheckLevels().first;
+    theProcess->SetEnergyMomentumCheckLevels(rellevel, abslevel);
+  }
+}
+
+
+void G4HadronicProcessStore::SetProcessRelLevel(G4double rellevel)
+{
+  G4cout << " Setting relative energy/momentum test level to " << rellevel << G4endl;
+  G4double abslevel = 0.0;
+  G4HadronicProcess* theProcess = 0;
+  for (G4int i = 0; i < G4int(process.size()); i++) {
+    theProcess = process[i];
+    abslevel = theProcess->GetEnergyMomentumCheckLevels().second;
+    theProcess->SetEnergyMomentumCheckLevels(rellevel, abslevel);
+  }
+}
