Index: trunk/source/processes/electromagnetic/lowenergy/src/G4AtomicDeexcitation.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4AtomicDeexcitation.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4AtomicDeexcitation.cc	(revision 1055)
@@ -26,5 +26,5 @@
 //
 // $Id: G4AtomicDeexcitation.cc,v 1.11 
-// GEANT4 tag $Name: geant4-09-02 $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Authors: Elena Guardincerri (Elena.Guardincerri@ge.infn.it)
@@ -111,5 +111,5 @@
     }
   
-  // Look this in a particular way: only one auger emitted! //
+  // Look this in a particular way: only one auger emitted! // ????
   while (provShellId > -2); 
   
@@ -297,6 +297,6 @@
 //  	      G4cout << "G4AtomicDeexcitation warning: No Auger transition found" <<  G4endl;
 // 	      G4cout << "Absorbed enrgy deposited locally" << G4endl;
- 	      return 0;
-// 	      //  G4Exception("G4AtomicDeexcitation: No Auger transition found");
+ 	      G4Exception("G4AtomicDeexcitation: No Auger transition found");
+	      return 0;
  	    }
 	}
@@ -394,6 +394,8 @@
       // G4int augerOriginatingShellId = 0;
       
-      G4int numberOfPossibleAuger = 
-	  (anAugerTransition->AugerTransitionProbabilities(transitionRandomShellId))->size();
+      G4int numberOfPossibleAuger = 0;
+      numberOfPossibleAuger = anAugerTransition->AugerTransitionProbabilities(transitionRandomShellId)->size();
+
+
       G4bool foundFlag = false;
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4AugerData.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4AugerData.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4AugerData.cc	(revision 1055)
@@ -486,5 +486,5 @@
       //	G4cout << "G4AugerData for Element no. " << element << " are loaded" << G4endl;
       // G4cout << "G4AugerData for Element no. " << element << " are loaded" << G4endl;
-      G4cout << "AugerTransitionTable complete"<< G4endl;
+      //G4cout << "AugerTransitionTable complete"<< G4endl;
     }
 }
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4BremsstrahlungCrossSectionHandler.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4BremsstrahlungCrossSectionHandler.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4BremsstrahlungCrossSectionHandler.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4BremsstrahlungCrossSectionHandler.cc,v 1.9 2006/06/29 19:38:42 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4BremsstrahlungCrossSectionHandler.cc,v 1.10 2009/03/03 11:19:18 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // -------------------------------------------------------------------
@@ -41,4 +41,6 @@
 // 10.10.2001 MGP Revision to improve code quality and consistency with design
 // 21.01.2003 VI  cut per region
+// 03.03.2009 LP  Added public method to make a easier migration of
+//                G4LowEnergyBremsstrahlung to G4LivermoreBremsstrahlungModel
 //
 // -------------------------------------------------------------------
@@ -55,4 +57,6 @@
 #include "G4ProductionCutsTable.hh"
 
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
+
 G4BremsstrahlungCrossSectionHandler::G4BremsstrahlungCrossSectionHandler(const G4VEnergySpectrum* spec,
 									 G4VDataSetAlgorithm* )
@@ -62,4 +66,5 @@
 }
 
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
 G4BremsstrahlungCrossSectionHandler::~G4BremsstrahlungCrossSectionHandler()
@@ -68,4 +73,5 @@
 }
 
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
 std::vector<G4VEMDataSet*>*
@@ -127,2 +133,18 @@
   return set;
 }
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
+
+G4double G4BremsstrahlungCrossSectionHandler::GetCrossSectionAboveThresholdForElement(G4double energy,
+										      G4double cutEnergy,
+										      G4int Z)
+{
+  G4double value = 0.;
+  if(energy > cutEnergy)
+    {
+      G4double elemCs = FindValue(Z, energy);
+      value  = theBR->Probability(Z,cutEnergy, energy, energy);	  
+      value *= elemCs;
+    }
+  return value;
+}
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeDecrease.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeDecrease.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeDecrease.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionChargeDecrease.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionChargeDecrease.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionChargeDecrease.hh"
@@ -79,4 +79,13 @@
   }
 
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionChargeDecrease is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeIncrease.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeIncrease.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionChargeIncrease.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionChargeIncrease.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionChargeIncrease.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionChargeIncrease.hh"
@@ -79,4 +79,13 @@
   }
 
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionChargeIncrease is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticChampion.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticChampion.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticChampion.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionElasticChampion.cc,v 1.4 2008/12/05 11:58:16 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionElasticChampion.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // -------------------------------------------------------------------
 
@@ -57,4 +57,14 @@
     G4Exception("G4CrossSectionElasticChampion constructor: electron is not defined");
   }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionElasticChampion is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordHE.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordHE.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordHE.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionElasticScreenedRutherfordHE.cc,v 1.2 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionElasticScreenedRutherfordHE.cc,v 1.3 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionElasticScreenedRutherfordHE.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 200. * eV; 
   highEnergyLimit = 10. * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionElasticScreenedRutherfordHE is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordLE.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordLE.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionElasticScreenedRutherfordLE.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionElasticScreenedRutherfordLE.cc,v 1.2 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionElasticScreenedRutherfordLE.cc,v 1.3 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionElasticScreenedRutherfordLE.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 0 * eV; 
   highEnergyLimit = 200 * eV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionElasticScreenedRutherfordLE is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationBorn.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationBorn.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationBorn.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionExcitationBorn.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionExcitationBorn.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionExcitationBorn.hh"
@@ -36,4 +36,14 @@
   highEnergyLimit = 10 * MeV;
   table = 0;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionExcitationBorn is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationEmfietzoglou.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationEmfietzoglou.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationEmfietzoglou.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionExcitationEmfietzoglou.cc,v 1.5 2008/12/05 11:58:16 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionExcitationEmfietzoglou.cc,v 1.6 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionExcitationEmfietzoglou.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 8.23 * eV;
   highEnergyLimit = 10. * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionExcitationEmfietzoglou is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreen.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreen.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreen.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionExcitationMillerGreen.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionExcitationMillerGreen.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionExcitationMillerGreen.hh"
@@ -92,4 +92,14 @@
   }
     
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionExcitationMillerGreen is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreenPartial.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreenPartial.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionExcitationMillerGreenPartial.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionExcitationMillerGreenPartial.cc,v 1.2 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionExcitationMillerGreenPartial.cc,v 1.3 2009/01/20 07:40:53 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionExcitationMillerGreenPartial.hh"
@@ -115,4 +115,8 @@
   tCorrected = k * kineticEnergyCorrection[particleTypeIndex];
 
+  // SI - added protection 
+  if (tCorrected < waterExcitation.ExcitationEnergy(excitationLevel)) return 0;
+  //
+  
   G4int z = 10;
 
@@ -158,5 +162,7 @@
   // ELECTRON CORRECTION
  
-  if ( particle == instance->GetIon("alpha++"))
+  if ( particle == instance->GetIon("alpha++")||
+       particle == G4Proton::ProtonDefinition())
+       
   {  while (i > 0)
      {
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBorn.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBorn.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBorn.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionIonisationBorn.cc,v 1.4 2008/08/20 14:51:48 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionIonisationBorn.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionIonisationBorn.hh"
@@ -82,4 +82,14 @@
     G4Exception("G4CrossSectionIonisationBorn Constructor: proton is not defined");
   }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionIonisationBorn is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBornPartial.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBornPartial.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationBornPartial.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionIonisationBornPartial.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionIonisationBornPartial.cc,v 1.5 2009/01/20 07:40:53 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionIonisationBornPartial.hh"
@@ -33,5 +33,5 @@
 G4CrossSectionIonisationBornPartial::G4CrossSectionIonisationBornPartial()
 {
-  lowEnergyLimitDefault = 25 * eV;
+  lowEnergyLimitDefault = 12.61 * eV;
   highEnergyLimitDefault = 30 * keV;
 
@@ -52,5 +52,5 @@
     tableFile[electron] = fileElectron;
 
-    lowEnergyLimit[electron] = 25. * eV;
+    lowEnergyLimit[electron] = 12.61 * eV;
     highEnergyLimit[electron] = 30. * keV;
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationRudd.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationRudd.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4CrossSectionIonisationRudd.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4CrossSectionIonisationRudd.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4CrossSectionIonisationRudd.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4CrossSectionIonisationRudd.hh"
@@ -154,4 +154,15 @@
     G4Exception("G4CrossSectionIonisationRudd Constructor: helium is not defined");
   }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4CrossSectionIonisationRudd is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
+
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeDecrease.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeDecrease.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeDecrease.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateChargeDecrease.cc,v 1.3 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateChargeDecrease.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateChargeDecrease.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 1 * keV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateChargeDecrease is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
@@ -70,7 +80,11 @@
   }
   
+  //SI - Added protection against total energy deposit
+  product.DoNotDepositEnergy();
+  //
   product.KillPrimaryParticle();
+
   product.AddEnergyDeposit(waterBindingEnergy);
-  
+
   G4DynamicParticle* aSecondary = new G4DynamicParticle(OutgoingParticleDefinition(definition, finalStateIndex), 
 							track.GetDynamicParticle()->GetMomentumDirection(), 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeIncrease.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeIncrease.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateChargeIncrease.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateChargeIncrease.cc,v 1.3 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateChargeIncrease.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateChargeIncrease.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 1 * keV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class 4FinalStateChargeIncrease is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
@@ -48,6 +58,10 @@
   product.Clear();
 
+  //SI - Added protection against total energy deposit
+  product.DoNotDepositEnergy();
+  //
   product.KillPrimaryParticle();
   product.AddEnergyDeposit(0.);
+
   G4ParticleDefinition* definition = track.GetDefinition();
  
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticBrennerZaider.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticBrennerZaider.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticBrennerZaider.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateElasticBrennerZaider.cc,v 1.8 2008/12/05 11:58:16 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateElasticBrennerZaider.cc,v 1.9 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateElasticBrennerZaider.hh"
@@ -64,4 +64,14 @@
   gamma100_200Coeff.push_back(-2.96264E-5);
   gamma100_200Coeff.push_back(-1.20655E-7);
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateElasticBrennerZaider is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticChampion.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticChampion.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticChampion.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateElasticChampion.cc,v 1.7 2008/12/10 18:25:28 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateElasticChampion.cc,v 1.8 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // -------------------------------------------------------------------
 
@@ -84,4 +84,14 @@
       G4Exception("G4FinalStateElastiChampion : constructor: electron is not defined");
   }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateElastiChampion is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticScreenedRutherford.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticScreenedRutherford.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateElasticScreenedRutherford.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateElasticScreenedRutherford.cc,v 1.4 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateElasticScreenedRutherford.cc,v 1.5 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateElasticScreenedRutherford.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 200 * eV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateElasticScreenedRutherford is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationBorn.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationBorn.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationBorn.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateExcitationBorn.cc,v 1.3 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateExcitationBorn.cc,v 1.4 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateExcitationBorn.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 500 * keV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateExcitationBorn is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationEmfietzoglou.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationEmfietzoglou.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationEmfietzoglou.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateExcitationEmfietzoglou.cc,v 1.5 2008/12/05 11:58:16 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateExcitationEmfietzoglou.cc,v 1.6 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateExcitationEmfietzoglou.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 8.23 * eV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateExcitationEmfietzoglou is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationMillerGreen.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationMillerGreen.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateExcitationMillerGreen.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateExcitationMillerGreen.cc,v 1.3 2008/07/14 20:47:34 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateExcitationMillerGreen.cc,v 1.4 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateExcitationMillerGreen.hh"
@@ -35,4 +35,14 @@
   lowEnergyLimit = 10 * eV;
   highEnergyLimit = 10 * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateExcitationMillerGreen is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationBorn.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationBorn.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationBorn.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateIonisationBorn.cc,v 1.16 2008/12/06 13:47:12 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateIonisationBorn.cc,v 1.17 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateIonisationBorn.hh"
@@ -124,4 +124,14 @@
     G4Exception("G4FinalStateIonisationBorn Constructor: proton is not defined");
   }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateIonisationBorn is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationRudd.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationRudd.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateIonisationRudd.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4FinalStateIonisationRudd.cc,v 1.8 2008/08/20 14:51:48 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateIonisationRudd.cc,v 1.9 2009/05/02 15:07:47 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 
 #include "G4FinalStateIonisationRudd.hh"
@@ -70,4 +70,14 @@
   lowEnergyLimit[helium] = 1. * keV;
   highEnergyLimit[helium] = 10. * MeV;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4FinalStateIonisationRudd is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateProduct.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateProduct.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4FinalStateProduct.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4FinalStateProduct.cc,v 1.5 2007/11/09 20:11:04 pia Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4FinalStateProduct.cc,v 1.6 2009/01/20 07:50:28 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // Contact Author: Maria Grazia Pia (Maria.Grazia.Pia@cern.ch)
@@ -53,5 +53,5 @@
 #include "G4ThreeVector.hh"
 
-G4FinalStateProduct::G4FinalStateProduct() : killStatus(false), isModified(false), localEnergyDeposit(0.), modifiedEnergy(0)
+G4FinalStateProduct::G4FinalStateProduct() : killStatus(false), doNotDepositStatus(false), isModified(false), localEnergyDeposit(0.), modifiedEnergy(0)
 {
   // empty
@@ -68,4 +68,5 @@
   // Reset object status
   killStatus = false;
+  doNotDepositStatus = false;
   isModified = false;
   localEnergyDeposit = 0.;
@@ -98,6 +99,12 @@
 }
   
+void G4FinalStateProduct::DoNotDepositEnergy()
+{
+  doNotDepositStatus = true;
+}
+
 void G4FinalStateProduct::KillPrimaryParticle()
 {
+  
   // ---- MGP ---- To be added: Handle local energy deposit here
   killStatus = true;
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4IonParametrisedLossModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4IonParametrisedLossModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4IonParametrisedLossModel.cc	(revision 1055)
@@ -37,6 +37,11 @@
 // First implementation: 10. 11. 2008
 //
-// Modifications:
-//
+// Modifications: 03. 02. 2009 - Bug fix iterators (AL)
+//                11. 03. 2009 - Introduced new table handler (G4IonDEDXHandler)
+//                               and modified method to add/remove tables
+//                               (tables are now built in initialisation phase),
+//                               Minor bug fix in ComputeDEDXPerVolume (AL)
+//                11. 05. 2009 - Introduced scaling algorithm for heavier ions:
+//                               G4IonDEDXScalingICRU73 (AL)
 //
 // Class description:
@@ -55,6 +60,11 @@
 #include "G4MaterialStoppingICRU73.hh"
 #include "G4SimpleMaterialStoppingICRU73.hh"
+#include "G4IronStoppingICRU73.hh"
+#include "G4VIonDEDXTable.hh"
+#include "G4VIonDEDXScalingAlgorithm.hh"
+#include "G4IonDEDXScalingICRU73.hh"
 #include "G4BraggIonModel.hh"
 #include "G4BetheBlochModel.hh"
+#include "G4ProductionCutsTable.hh"
 #include "G4ParticleChangeForLoss.hh"
 #include "G4LossTableManager.hh"
@@ -63,4 +73,8 @@
 #include "Randomize.hh"
 
+//#define PRINT_TABLE_BUILT
+
+
+// #########################################################################
 
 G4IonParametrisedLossModel::G4IonParametrisedLossModel(
@@ -89,7 +103,44 @@
   betheBlochModel = new G4BetheBlochModel();
 
-  // By default ICRU 73 stopping power tables are loaded
-  AddDEDXTable<G4SimpleMaterialStoppingICRU73>();
-  AddDEDXTable<G4MaterialStoppingICRU73>();
+  // By default ICRU 73 stopping power tables are loaded:
+
+  // Ions with Z above between 19 and 21: Ar-40 data is used as basis 
+  // for stopping power scaling
+  G4int ionZMin = 19;
+  G4int ionZMax = 21;
+  G4int refIonZ = 18;
+  G4int refIonA = 40;
+
+  AddDEDXTable("ICRU73-elemmat",
+           new G4SimpleMaterialStoppingICRU73,
+	   new G4IonDEDXScalingICRU73(ionZMin, ionZMax, refIonZ, refIonA));
+
+  // Ions with Z above 21: Fe-56 data is used as basis for stopping power 
+  // scaling
+  ionZMin = 22;
+  ionZMax = 102;
+  refIonZ = 26;
+  refIonA = 56;
+
+  AddDEDXTable("ICRU73-ironions",
+               new G4IronStoppingICRU73,
+	       new G4IonDEDXScalingICRU73(ionZMin, ionZMax, refIonZ, refIonA));
+
+  // Compound materials: Ar-40 data is used as basis for stopping power 
+  // scaling (except for iron ions)
+  ionZMin = 19;
+  ionZMax = 102;
+  refIonZ = 18;
+  refIonA = 40;
+
+  G4IonDEDXScalingICRU73* scaling =
+            new G4IonDEDXScalingICRU73(ionZMin, ionZMax, refIonZ, refIonA);
+
+  G4int ironIonAtomicNumber = 26;
+  scaling -> AddException(ironIonAtomicNumber);
+
+  AddDEDXTable("ICRU73-compmat",
+               new G4MaterialStoppingICRU73,
+	       scaling);
 
   // The boundaries for the range tables are set
@@ -107,5 +158,5 @@
   dedxCacheMaterial = 0;
   dedxCacheEnergyCut = 0;
-  dedxCacheIter = lossTableList.begin();
+  dedxCacheIter = lossTableList.end();
   dedxCacheTransitionEnergy = 0.0;  
   dedxCacheTransitionFactor = 0.0;
@@ -113,4 +164,5 @@
 }
 
+// #########################################################################
 
 G4IonParametrisedLossModel::~G4IonParametrisedLossModel() {
@@ -142,4 +194,5 @@
 }
 
+// #########################################################################
 
 G4double G4IonParametrisedLossModel::MinEnergyCut(
@@ -151,4 +204,5 @@
 }
 
+// #########################################################################
 
 void G4IonParametrisedLossModel::Initialise(
@@ -166,5 +220,5 @@
   dedxCacheMaterial = 0;
   dedxCacheEnergyCut = 0;
-  dedxCacheIter = lossTableList.begin();
+  dedxCacheIter = lossTableList.end();
   dedxCacheTransitionEnergy = 0.0;  
   dedxCacheTransitionFactor = 0.0;
@@ -196,4 +250,52 @@
   cutEnergies.clear();
   for(size_t i = 0; i < size; i++) cutEnergies.push_back(cuts[i]);
+
+  // All dE/dx vectors are built
+  const G4ProductionCutsTable* coupleTable=
+                     G4ProductionCutsTable::GetProductionCutsTable();
+  size_t nmbCouples = coupleTable -> GetTableSize();
+
+#ifdef PRINT_TABLE_BUILT
+    G4cout << "G4IonParametrisedLossModel::Initialise():"
+           << " Building dE/dx vectors:"
+           << G4endl;         
+#endif
+
+  for (size_t i = 0; i < nmbCouples; i++) {
+
+    const G4MaterialCutsCouple* couple = 
+                                     coupleTable -> GetMaterialCutsCouple(i);
+
+    const G4Material* material = couple -> GetMaterial();
+    //    G4ProductionCuts* productionCuts = couple -> GetProductionCuts();
+
+    for(G4int atomicNumberIon = 3; atomicNumberIon < 102; atomicNumberIon++) {
+
+       LossTableList::iterator iter = lossTableList.begin();
+       LossTableList::iterator iter_end = lossTableList.end();
+
+       for(;iter != iter_end; iter++) { 
+
+          if(*iter == 0) {
+              G4cout << "G4IonParametrisedLossModel::Initialise():"
+                     << " Skipping illegal table."
+                     << G4endl;         
+          }
+
+          G4bool isApplicable = 
+                    (*iter) -> BuildDEDXTable(atomicNumberIon, material);
+          if(isApplicable) {
+
+#ifdef PRINT_TABLE_BUILT
+             G4cout << "  Atomic Number Ion = " << atomicNumberIon
+                    << ", Material = " << material -> GetName()
+                    << ", Table = " << (*iter) -> GetName()
+                    << G4endl;      
+#endif
+             break; 
+	  }
+       }
+    }
+  }
 
   // The particle change object is cast to G4ParticleChangeForLoss
@@ -221,4 +323,5 @@
 }
 
+// #########################################################################
 
 G4double G4IonParametrisedLossModel::ComputeCrossSectionPerAtom(
@@ -285,4 +388,5 @@
 }
 
+// #########################################################################
 
 G4double G4IonParametrisedLossModel::CrossSectionPerVolume(
@@ -303,4 +407,5 @@
 }
 
+// #########################################################################
 
 G4double G4IonParametrisedLossModel::ComputeDEDXPerVolume(
@@ -342,5 +447,5 @@
   LossTableList::iterator iter = dedxCacheIter;
 
-  if(iter != lossTableList.begin()) {
+  if(iter != lossTableList.end()) {
 
      G4double transitionEnergy = dedxCacheTransitionEnergy;
@@ -369,17 +474,18 @@
         if(scaledTransitionEnergy >= lowEnergyLimit) {
 
-           G4double factor = 1.0 + dedxCacheTransitionFactor / 
-                                                           kineticEnergy;
-
            dEdx = betheBlochModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
      				      scaledKineticEnergy, cutEnergy);
+        
+           dEdx *= chargeSquare;
+
+           dEdx += corrections -> ComputeIonCorrections(particle, 
+                                                 material, kineticEnergy);
+
+           G4double factor = 1.0 + dedxCacheTransitionFactor / 
+                                                           kineticEnergy;
+
 	   dEdx *= factor;
-
-        }
-        dEdx *= chargeSquare;
-
-        dEdx += corrections -> ComputeIonCorrections(particle, 
-                                                 material, kineticEnergy);
+	}
      }
   }
@@ -446,36 +552,8 @@
   if (dEdx < 0.0) dEdx = 0.0;
 
-#ifdef PRINT_DEBUG
-
-  G4cout << "########################################################"
-         << G4endl
-         << "# G4IonParametrisedLossModel::ComputeDEDXPerVolume"
-         << G4endl
-         << "# Material =" << material -> GetName() 
-         << G4endl
-         << "# Particle = " << particle -> GetParticleName()            
-         << G4endl;
-         << "# Cut energy (MeV) = " << cutEnergy/MeV  
-         << G4endl;
-  
-  G4cout << "#"
-         << std::setw(13) << std::right << "E(MeV)"
-         << std::setw(14) << "dE/dx(keV/um)"
-         << std::setw(14) << "d:dE/dx(keV/um)"
-         << std::setw(14) << "(d:dE/dx)/dE/dx"
-         << G4endl
-         << "# ------------------------------------------------------"
-         << G4endl;
-
-  G4cout << std::setw(14) << std::right << kineticEnergy / MeV
-         << std::setw(14) << (dEdx + dEdXDeltaRays) / keV * um
-         << std::setw(14) << dEdXDeltaRays / keV * um
-         << std::setw(14) << dEdXDeltaRays / (dEdx + dEdXDeltaRays) * 100.0
-         << G4endl;
-#endif
-
   return dEdx;
 }
 
+// #########################################################################
 
 void G4IonParametrisedLossModel::PrintDEDXTable(
@@ -508,5 +586,5 @@
          << std::setw(13) << std::right << "(MeV)"
          << std::setw(14) << "(MeV)"
-         << std::setw(14) << "(MeV/mm)"
+         << std::setw(14) << "(MeV/cm)"
          << std::setw(14) << "(MeV*cm2/mg)"
          << G4endl
@@ -535,5 +613,5 @@
       G4cout << std::setw(14) << std::right << energy / MeV
              << std::setw(14) << energy / atomicMassNumber / MeV
-             << std::setw(14) << dedx / MeV * mm
+             << std::setw(14) << dedx / MeV * cm
              << std::setw(14) << dedx / materialDensity / (MeV*cm2/(0.001*g)) 
              << G4endl;
@@ -541,4 +619,29 @@
 }
 
+// #########################################################################
+
+void G4IonParametrisedLossModel::PrintDEDXTableHandlers(
+                   const G4ParticleDefinition* particle,  // Projectile (ion) 
+                   const G4Material* material,  // Absorber material
+                   G4double lowerBoundary,      // Minimum energy per nucleon
+                   G4double upperBoundary,      // Maximum energy per nucleon
+                   G4int nmbBins,               // Number of bins
+                   G4bool logScaleEnergy) {     // Logarithmic scaling of energy
+
+  LossTableList::iterator iter = lossTableList.begin();
+  LossTableList::iterator iter_end = lossTableList.end();
+
+  for(;iter != iter_end; iter++) { 
+      G4bool isApplicable =  (*iter) -> IsApplicable(particle, material);
+      if(isApplicable) {  
+	(*iter) -> PrintDEDXTable(particle, material,
+                                  lowerBoundary, upperBoundary, 
+                                  nmbBins,logScaleEnergy); 
+        break;
+      } 
+  }
+}
+
+// #########################################################################
 
 void G4IonParametrisedLossModel::SampleSecondaries(
@@ -641,4 +744,5 @@
 }
 
+// #########################################################################
 
 void G4IonParametrisedLossModel::UpdateDEDXCache(
@@ -672,5 +776,5 @@
 
      // If any table is applicable, the transition factor is computed:
-     if(iter != lossTableList.begin()) {
+     if(iter != lossTableList.end()) {
 
         // Retrieving the transition energy from the parameterisation table
@@ -740,4 +844,5 @@
 }
 
+// #########################################################################
 
 void G4IonParametrisedLossModel::CorrectionsAlongStep(
@@ -767,4 +872,6 @@
   G4double kineticEnergy = dynamicParticle -> GetKineticEnergy();
 
+  if(kineticEnergy == eloss) { return; }
+
   G4double cutEnergy = DBL_MAX;
   size_t cutIndex = couple -> GetIndex();
@@ -777,5 +884,5 @@
   // If parameterization for ions is available the electronic energy loss
   // is overwritten
-  if(iter != lossTableList.begin()) {
+  if(iter != lossTableList.end()) {
 
      // The energy loss is calculated using the ComputeDEDXPerVolume function
@@ -832,5 +939,4 @@
 
      }
-
   }
 
@@ -854,5 +960,5 @@
   G4double transitionEnergy = dedxCacheTransitionEnergy;
 
-  if(iter != lossTableList.begin() && transitionEnergy < kineticEnergy) {
+  if(iter != lossTableList.end() && transitionEnergy < kineticEnergy) {
      chargeSquareRatio *= corrections -> EffectiveChargeCorrection(particle,
                                                                 material,
@@ -862,5 +968,5 @@
      eloss *= chargeSquareRatioCorr;
   }
-  else if (iter == lossTableList.begin()) {
+  else if (iter == lossTableList.end()) {
 
      chargeSquareRatio *= corrections -> EffectiveChargeCorrection(particle,
@@ -875,5 +981,5 @@
   // overwrite the energy loss (i.e. when the effective charge approach is
   // used)
-  if(iter == lossTableList.begin()) {
+  if(iter == lossTableList.end()) {
 
      G4double scaledKineticEnergy = kineticEnergy * dedxCacheGenIonMassRatio;
@@ -906,7 +1012,7 @@
      particleChangeLoss -> ProposeNonIonizingEnergyDeposit(nloss);
   } 
-
-}
-
+}
+
+// #########################################################################
 
 void G4IonParametrisedLossModel::BuildRangeVector(
@@ -1018,9 +1124,58 @@
 
   IonMatCouple ionMatCouple = std::make_pair(particle, material); 
- 
+
   E[ionMatCouple] = energyRangeVector;
   r[ionMatCouple] = rangeEnergyVector;
 }
 
+// #########################################################################
+
+G4double G4IonParametrisedLossModel::GetRange(
+                  const G4ParticleDefinition* particle, // Projectile
+		  const G4Material* material,           // Target Material
+	          G4double kineticEnergy) {
+
+  G4double range = 0.0;
+
+  IonMatCouple couple = std::make_pair(particle, material);
+
+  EnergyRangeTable::iterator iter = E.find(couple);
+
+  if(iter == E.end()) {
+     G4cerr << "G4IonParametrisedLossModel::GetRange() No range vector found."
+            << G4endl;
+
+     G4cout << "   Ion-material pair: " << particle ->GetParticleName()
+            << "  " << material -> GetName()
+            << G4endl
+            << "   Available couples:"
+            << G4endl;
+
+     EnergyRangeTable::iterator iter_beg = E.begin();
+     EnergyRangeTable::iterator iter_end = E.end();
+
+     for(;iter_beg != iter_end; iter_beg++) {
+         IonMatCouple key = (*iter_beg).first; 
+ 
+         G4cout << "           " << (key.first) -> GetParticleName()
+                << "  " << (key.second) -> GetName()
+                << G4endl;
+     }
+  }
+  else {
+     G4PhysicsVector* energyRange = (*iter).second;
+
+     if(energyRange != 0) {
+        G4bool b;
+
+        // Computing range for kinetic energy:
+        range = energyRange -> GetValue(kineticEnergy, b);
+     }
+  }
+
+  return range;
+}
+
+// #########################################################################
 
 G4double G4IonParametrisedLossModel::ComputeLossForStep(
@@ -1061,18 +1216,74 @@
         if(loss < 0.0) loss = 0.0;
      }
-
-#ifdef PRINT_DEBUG
-     G4cout << "G4IonParametrisedLossModel::ComputeLossForStep() E = " 
-            << kineticEnergy / MeV << " MeV * "
-            << value.energyScaling << " = "
-            << kineticEnergy * value.energyScaling / MeV 
-            << " MeV, dE/dx = " << dedx / MeV * cm << " MeV/cm = "
-            << dedx/factor/MeV*cm << " * " << factor << " MeV/cm; index = "
-            << value.dEdxIndex << ", material = " << material -> GetName()
-            << G4endl;
-#endif
-
   }
 
   return loss;
 }
+
+// #########################################################################
+
+G4bool G4IonParametrisedLossModel::AddDEDXTable(
+                                const G4String& name,
+                                G4VIonDEDXTable* table, 
+                                G4VIonDEDXScalingAlgorithm* algorithm) {
+
+  if(table == 0) {
+     G4cerr << "G4IonParametrisedLossModel::AddDEDXTable() Cannot "
+            << " add table: Invalid pointer."
+            << G4endl;      
+
+     return false;
+  }
+
+  // Checking uniqueness of name
+  LossTableList::iterator iter = lossTableList.begin();
+  LossTableList::iterator iter_end = lossTableList.end();
+  
+  for(;iter != iter_end; iter++) {
+     G4String tableName = (*iter) -> GetName();
+
+     if(tableName == name) { 
+        G4cerr << "G4IonParametrisedLossModel::AddDEDXTable() Cannot "
+               << " add table: Name already exists."      
+               << G4endl;
+
+        return false;
+     }
+  }
+
+  G4VIonDEDXScalingAlgorithm* scalingAlgorithm = algorithm;
+  if(scalingAlgorithm == 0) 
+     scalingAlgorithm = new G4VIonDEDXScalingAlgorithm; 
+ 
+  G4IonDEDXHandler* handler = 
+                      new G4IonDEDXHandler(table, scalingAlgorithm, name); 
+
+  lossTableList.push_front(handler);
+
+  return true;
+}
+
+// #########################################################################
+
+G4bool G4IonParametrisedLossModel::RemoveDEDXTable(
+				 const G4String& name) {
+
+  LossTableList::iterator iter = lossTableList.begin();
+  LossTableList::iterator iter_end = lossTableList.end();
+  
+  for(;iter != iter_end; iter++) {
+     G4String tableName = (*iter) -> GetName();
+
+     if(tableName == name) { 
+        delete (*iter);
+
+        lossTableList.erase(iter);
+        return true;
+     }
+  }
+
+  return false;
+}
+
+// #########################################################################
+
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreComptonModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreComptonModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreComptonModel.cc	(revision 1055)
@@ -24,7 +24,19 @@
 // ********************************************************************
 //
-// $Id: G4LivermoreComptonModel.cc,v 1.1 2008/10/30 14:17:46 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
-//
+// $Id: G4LivermoreComptonModel.cc,v 1.6 2009/04/18 18:29:34 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
+//
+//
+// Author: Sebastien Inserti
+//         30 October 2008
+//
+// History:
+// --------
+// 18 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - remove GetMeanFreePath method and table
+//                  - added protection against numerical problem in energy sampling 
+//                  - use G4ElementSelector
 
 #include "G4LivermoreComptonModel.hh"
@@ -37,13 +49,14 @@
 
 G4LivermoreComptonModel::G4LivermoreComptonModel(const G4ParticleDefinition*,
-                                             const G4String& nam)
-:G4VEmModel(nam),isInitialised(false)
+						 const G4String& nam)
+  :G4VEmModel(nam),isInitialised(false),meanFreePathTable(0),
+   scatterFunctionData(0),crossSectionHandler(0)
 {
-  lowEnergyLimit = 250 * eV; // SI - Could be 10 eV ?
+  lowEnergyLimit = 250 * eV; 
   highEnergyLimit = 100 * GeV;
-  SetLowEnergyLimit(lowEnergyLimit);
+  //  SetLowEnergyLimit(lowEnergyLimit);
   SetHighEnergyLimit(highEnergyLimit);
 
-  verboseLevel= 0;
+  verboseLevel=0 ;
   // Verbosity scale:
   // 0 = nothing 
@@ -52,11 +65,12 @@
   // 3 = calculation of cross sections, file openings, sampling of atoms
   // 4 = entering in methods
-  
-  G4cout << "Livermore Compton model is constructed " << G4endl
-         << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
-         << highEnergyLimit / GeV << " GeV"
-         << G4endl;
- 
+
+  if(  verboseLevel>0 ) { 
+    G4cout << "Livermore Compton model is constructed " << G4endl
+	   << "Energy range: "
+	   << lowEnergyLimit / eV << " eV - "
+	   << highEnergyLimit / GeV << " GeV"
+	   << G4endl;
+  }
 }
 
@@ -65,7 +79,6 @@
 G4LivermoreComptonModel::~G4LivermoreComptonModel()
 {  
-  delete meanFreePathTable;
-  delete crossSectionHandler;
-  delete scatterFunctionData;
+  if (crossSectionHandler) delete crossSectionHandler;
+  if (scatterFunctionData) delete scatterFunctionData;
 }
 
@@ -73,27 +86,15 @@
 
 void G4LivermoreComptonModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+					 const G4DataVector& cuts)
 {
   if (verboseLevel > 3)
     G4cout << "Calling G4LivermoreComptonModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
-  
-  if (LowEnergyLimit() < lowEnergyLimit)
-    {
-      G4cout << "G4LivermoreComptonModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
-      SetLowEnergyLimit(lowEnergyLimit);
-    }
-
-  if (HighEnergyLimit() > highEnergyLimit)
-    {
-      G4cout << "G4LivermoreComptonModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
-      SetHighEnergyLimit(highEnergyLimit);
-    }
-
+  if (crossSectionHandler)
+  {
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
+  }
+  
   // Reading of data files - all materials are read
   
@@ -113,25 +114,19 @@
   shellData.LoadData(file);
 
-  meanFreePathTable = 0;
-  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
-  
   if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for Livermore Compton model" << G4endl;
 
-  G4cout << "Livermore Compton model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
-
-  //
-  
+  InitialiseElementSelectors(particle,cuts);
+
+  if(  verboseLevel>0 ) { 
+    G4cout << "Livermore Compton model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / eV << " eV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
+  //  
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-    
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -148,5 +143,7 @@
     G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermoreComptonModel" << G4endl;
 
-  G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit) return 0.0;
+    
+  G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);  
   return cs;
 }
@@ -155,9 +152,9 @@
 
 void G4LivermoreComptonModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
-					      const G4MaterialCutsCouple* couple,
-					      const G4DynamicParticle* aDynamicGamma,
-					      G4double,
-					      G4double)
+						const G4MaterialCutsCouple* couple,
+						const G4DynamicParticle* aDynamicGamma,
+						G4double, G4double)
 {
+
   // The scattered gamma energy is sampled according to Klein - Nishina formula.
   // then accepted or rejected depending on the Scattering Function multiplied
@@ -174,17 +171,20 @@
   // (Nucl Phys 20(1960),15).
 
-  if (verboseLevel > 3)
-    G4cout << "Calling SampleSecondaries() of G4LivermoreComptonModel" << G4endl;
-
   G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
-  
-  if (photonEnergy0 <= lowEnergyLimit)
-  {
+
+  if (verboseLevel > 3) {
+    G4cout << "G4LivermoreComptonModel::SampleSecondaries() E(MeV)= " 
+	   << photonEnergy0/MeV << " in " << couple->GetMaterial()->GetName() 
+	   << G4endl;
+  }
+  
+  // low-energy gamma is absorpted by this process
+  if (photonEnergy0 <= lowEnergyLimit) 
+    {
       fParticleChange->ProposeTrackStatus(fStopAndKill);
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
-      // SI - IS THE FOLLOWING RETURN NECESSARY ?
       return ;
-  }
+    }
 
   G4double e0m = photonEnergy0 / electron_mass_c2 ;
@@ -192,5 +192,8 @@
 
   // Select randomly one element in the current material
-  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  //  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+  const G4Element* elm = SelectRandomAtom(couple,particle,photonEnergy0);
+  G4int Z = (G4int)elm->GetZ();
 
   G4double epsilon0 = 1. / (1. + 2. * e0m);
@@ -212,5 +215,6 @@
       if ( alpha1/(alpha1+alpha2) > G4UniformRand())
       {
-        epsilon = std::exp(-alpha1 * G4UniformRand());  // std::pow(epsilon0,G4UniformRand())
+	// std::pow(epsilon0,G4UniformRand())
+        epsilon = std::exp(-alpha1 * G4UniformRand());  
         epsilonSq = epsilon * epsilon;
       }
@@ -238,6 +242,6 @@
   // Doppler broadening -  Method based on:
   // Y. Namito, S. Ban and H. Hirayama, 
-  // "Implementation of the Doppler Broadening of a Compton-Scattered Photon Into the EGS4 Code" 
-  // NIM A 349, pp. 489-494, 1994
+  // "Implementation of the Doppler Broadening of a Compton-Scattered Photon 
+  // into the EGS4 Code", NIM A 349, pp. 489-494, 1994
   
   // Maximum number of sampling iterations
@@ -257,5 +261,6 @@
       eMax = photonEnergy0 - bindingE;
      
-      // Randomly sample bound electron momentum (memento: the data set is in Atomic Units)
+      // Randomly sample bound electron momentum 
+      // (memento: the data set is in Atomic Units)
       G4double pSample = profileData.RandomSelectMomentum(Z,shell);
       // Rescale from atomic units
@@ -299,15 +304,22 @@
 
   if (photonEnergy1 > 0.)
-  {
-    fParticleChange->SetProposedKineticEnergy(photonEnergy1) ;
-  }
+    {
+      fParticleChange->SetProposedKineticEnergy(photonEnergy1) ;
+    }
   else
-  {
-    fParticleChange->SetProposedKineticEnergy(0.) ;
-    fParticleChange->ProposeTrackStatus(fStopAndKill);
-  }
+    {
+      photonEnergy1 = 0.;
+      fParticleChange->SetProposedKineticEnergy(0.) ;
+      fParticleChange->ProposeTrackStatus(fStopAndKill);   
+    }
 
   // Kinematics of the scattered electron
   G4double eKineticEnergy = photonEnergy0 - photonEnergy1 - bindingE;
+
+  // protection against negative final energy: no e- is created
+  if(eKineticEnergy < 0.0) {
+    fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0 - photonEnergy1);
+    return;
+  }
   G4double eTotalEnergy = eKineticEnergy + electron_mass_c2;
 
@@ -329,28 +341,11 @@
   eDirection.rotateUz(photonDirection0);
 
-// SI - The range test has been removed wrt original G4LowEnergyCompton class
+  // SI - The range test has been removed wrt original G4LowEnergyCompton class
 
   fParticleChange->ProposeLocalEnergyDeposit(bindingE);
   
-  G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),eDirection,eKineticEnergy) ;
+  G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),
+						 eDirection,eKineticEnergy) ;
   fvect->push_back(dp);
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermoreComptonModel::GetMeanFreePath(const G4Track& track,
-					     G4double, // previousStepSize
-					     G4ForceCondition*)
-{
-  const G4DynamicParticle* photon = track.GetDynamicParticle();
-  G4double energy = photon->GetKineticEnergy();
-  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
-  size_t materialIndex = couple->GetIndex();
-
-  G4double meanFreePath;
-  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
-  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
-  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
-  return meanFreePath;
-}
-
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreGammaConversionModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreGammaConversionModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreGammaConversionModel.cc	(revision 1055)
@@ -24,7 +24,19 @@
 // ********************************************************************
 //
-// $Id: G4LivermoreGammaConversionModel.cc,v 1.1 2008/10/30 14:16:35 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
-//
+// $Id: G4LivermoreGammaConversionModel.cc,v 1.6 2009/05/02 09:14:43 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
+//
+//
+// Author: Sebastien Inserti
+//         30 October 2008
+//
+// History:
+// --------
+// 12 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - use CLHEP electron mass for low-enegry limit
+//                  - remove MeanFreePath method and table
+
 
 #include "G4LivermoreGammaConversionModel.hh"
@@ -37,16 +49,12 @@
 
 G4LivermoreGammaConversionModel::G4LivermoreGammaConversionModel(const G4ParticleDefinition*,
-                                             const G4String& nam)
-:G4VEmModel(nam),smallEnergy(2.*MeV),isInitialised(false)
-{
-  lowEnergyLimit = 1.022000 * MeV;
+								 const G4String& nam)
+  :G4VEmModel(nam),smallEnergy(2.*MeV),isInitialised(false),
+   crossSectionHandler(0),meanFreePathTable(0)
+{
+  lowEnergyLimit = 2.0*electron_mass_c2;
   highEnergyLimit = 100 * GeV;
-  
-  G4cout << "Livermore Gamma conversion is constructed " << G4endl
-         << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
-         << highEnergyLimit / GeV << " GeV"
-         << G4endl;
-	 
+  SetHighEnergyLimit(highEnergyLimit);
+  	 
   verboseLevel= 0;
   // Verbosity scale:
@@ -57,4 +65,11 @@
   // 4 = entering in methods
 
+  if(verboseLevel > 0) {
+    G4cout << "Livermore Gamma conversion is constructed " << G4endl
+	   << "Energy range: "
+	   << lowEnergyLimit / MeV << " MeV - "
+	   << highEnergyLimit / GeV << " GeV"
+	   << G4endl;
+  }
 }
 
@@ -63,44 +78,29 @@
 G4LivermoreGammaConversionModel::~G4LivermoreGammaConversionModel()
 {  
-  delete meanFreePathTable;
-  delete crossSectionHandler;
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-void G4LivermoreGammaConversionModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+  if (crossSectionHandler) delete crossSectionHandler;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
+
+void 
+G4LivermoreGammaConversionModel::Initialise(const G4ParticleDefinition*,
+					    const G4DataVector&)
 {
   if (verboseLevel > 3)
     G4cout << "Calling G4LivermoreGammaConversionModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
-  
-  if (LowEnergyLimit() < lowEnergyLimit)
+  if (crossSectionHandler)
   {
-    G4cout << "G4LivermoreGammaConversionModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
-    SetLowEnergyLimit(lowEnergyLimit);
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
   }
 
-  if (HighEnergyLimit() > highEnergyLimit)
-  {
-    G4cout << "G4LivermoreGammaConversionModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
-    SetHighEnergyLimit(highEnergyLimit);
-  }
-
   // Read data tables for all materials
   
   crossSectionHandler = new G4CrossSectionHandler();
-  crossSectionHandler->Initialise(0,1.0220*MeV,100.*GeV,400);
+  crossSectionHandler->Initialise(0,lowEnergyLimit,100.*GeV,400);
   G4String crossSectionFile = "pair/pp-cs-";
   crossSectionHandler->LoadData(crossSectionFile);
 
-  meanFreePathTable = 0;
-  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
-
   //
   
@@ -108,28 +108,30 @@
     G4cout << "Loaded cross section files for PenelopeGammaConversion" << G4endl;
 
-  G4cout << "Livermore Gamma Conversion model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / MeV << " MeV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if (verboseLevel > 0) { 
+    G4cout << "Livermore Gamma Conversion model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / MeV << " MeV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-  isInitialised = true;}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermoreGammaConversionModel::ComputeCrossSectionPerAtom(
-                                       const G4ParticleDefinition*,
-                                             G4double GammaEnergy,
-                                             G4double Z, G4double,
-                                             G4double, G4double)
-{
-  if (verboseLevel > 3)
-    G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermoreGammaConversionModel" << G4endl;
+  fParticleChange = GetParticleChangeForGamma();
+  isInitialised = true;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
+
+G4double 
+G4LivermoreGammaConversionModel::ComputeCrossSectionPerAtom(const G4ParticleDefinition*,
+							    G4double GammaEnergy,
+							    G4double Z, G4double,
+							    G4double, G4double)
+{
+  if (verboseLevel > 3) {
+    G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermoreGammaConversionModel" 
+	   << G4endl;
+  }
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit) return 0;
 
   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
@@ -173,14 +175,20 @@
     {
       // Select randomly one element in the current material
-      const G4Element* element = crossSectionHandler->SelectRandomElement(couple,photonEnergy);
+      //const G4Element* element = crossSectionHandler->SelectRandomElement(couple,photonEnergy);
+      const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+      const G4Element* element = SelectRandomAtom(couple,particle,photonEnergy);
 
       if (element == 0)
 	{
-	  G4cout << "G4LivermoreGammaConversionModel::SampleSecondaries - element = 0" << G4endl;
+	  G4cout << "G4LivermoreGammaConversionModel::SampleSecondaries - element = 0" 
+		 << G4endl;
+	  return;
 	}
       G4IonisParamElm* ionisation = element->GetIonisation();
-       if (ionisation == 0)
+      if (ionisation == 0)
 	{
-	  G4cout << "G4LivermoreGammaConversionModel::SampleSecondaries - ionisation = 0" << G4endl;
+	  G4cout << "G4LivermoreGammaConversionModel::SampleSecondaries - ionisation = 0" 
+		 << G4endl;
+	  return;
 	}
 
@@ -272,8 +280,7 @@
   // distribution with respect to the Z axis along the parent photon
   
-//  aParticleChange.SetNumberOfSecondaries(2) ; 
   G4double electronKineEnergy = std::max(0.,electronTotEnergy - electron_mass_c2) ;
   
-// SI - The range test has been removed wrt original G4LowEnergyGammaconversion class
+  // SI - The range test has been removed wrt original G4LowEnergyGammaconversion class
 
   G4ThreeVector electronDirection (dxEle, dyEle, dzEle);
@@ -287,5 +294,5 @@
   G4double positronKineEnergy = std::max(0.,positronTotEnergy - electron_mass_c2) ;
 
-// SI - The range test has been removed wrt original G4LowEnergyGammaconversion class
+  // SI - The range test has been removed wrt original G4LowEnergyGammaconversion class
 
   G4ThreeVector positronDirection (dxPos, dyPos, dzPos);
@@ -337,20 +344,2 @@
 } 
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermoreGammaConversionModel::GetMeanFreePath(const G4Track& track, 
-						     G4double, // previousStepSize
-						     G4ForceCondition*)
-{
-  const G4DynamicParticle* photon = track.GetDynamicParticle();
-  G4double energy = photon->GetKineticEnergy();
-  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
-  size_t materialIndex = couple->GetIndex();
-
-  G4double meanFreePath;
-  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
-  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
-  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
-  return meanFreePath;
-}
-
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePhotoElectricModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePhotoElectricModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePhotoElectricModel.cc	(revision 1055)
@@ -24,7 +24,20 @@
 // ********************************************************************
 //
-// $Id: G4LivermorePhotoElectricModel.cc,v 1.1 2008/10/30 14:16:35 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
-//
+// $Id: G4LivermorePhotoElectricModel.cc,v 1.7 2009/04/18 18:29:34 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
+//
+//
+// Author: Sebastien Inserti
+//         30 October 2008
+//
+// History:
+// --------
+// 15 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - remove GetMeanFreePath method and table
+//                  - simplify sampling of deexcitation by using cut in energy
+//                  - added protection against numerical problem in energy sampling 
+//                  - use G4ElementSelector
 
 #include "G4LivermorePhotoElectricModel.hh"
@@ -38,17 +51,14 @@
 G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel(const G4ParticleDefinition*,
                                              const G4String& nam)
-:G4VEmModel(nam),isInitialised(false)
-{
-  lowEnergyLimit = 250 * eV; // SI - Could be 10 eV ?
+:G4VEmModel(nam),isInitialised(false),meanFreePathTable(0),
+ crossSectionHandler(0),shellCrossSectionHandler(0),ElectronAngularGenerator(0)
+{
+  lowEnergyLimit = 250 * eV; 
   highEnergyLimit = 100 * GeV;
-  SetLowEnergyLimit(lowEnergyLimit);
+  //  SetLowEnergyLimit(lowEnergyLimit);
   SetHighEnergyLimit(highEnergyLimit);
-  
-  G4cout << "Livermore Compton is constructed " << G4endl
-         << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
-         << highEnergyLimit / GeV << " GeV"
-         << G4endl;
- 
+
+  ActivateAuger(false);
+   
   verboseLevel= 0;
   // Verbosity scale:
@@ -58,4 +68,11 @@
   // 3 = calculation of cross sections, file openings, sampling of atoms
   // 4 = entering in methods
+  if(verboseLevel>0) {
+    G4cout << "Livermore PhotoElectric is constructed " << G4endl
+	   << "Energy range: "
+	   << lowEnergyLimit / eV << " eV - "
+	   << highEnergyLimit / GeV << " GeV"
+	   << G4endl;
+  }
 }
 
@@ -64,38 +81,30 @@
 G4LivermorePhotoElectricModel::~G4LivermorePhotoElectricModel()
 {  
-  delete meanFreePathTable;
-  delete crossSectionHandler;
-  delete shellCrossSectionHandler;
-  delete ElectronAngularGenerator;
-}
-
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-void G4LivermorePhotoElectricModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+  if (crossSectionHandler) delete crossSectionHandler;
+  if (shellCrossSectionHandler) delete shellCrossSectionHandler;
+  if (ElectronAngularGenerator) delete ElectronAngularGenerator;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
+
+void 
+G4LivermorePhotoElectricModel::Initialise(const G4ParticleDefinition*,
+					  const G4DataVector&)
 {
   if (verboseLevel > 3)
     G4cout << "Calling G4LivermorePhotoElectricModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
-  
-  if (LowEnergyLimit() < lowEnergyLimit)
+  if (crossSectionHandler)
   {
-      G4cout << "G4LivermorePhotoElectricModel: low energy limit increased from " <<
-        LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << 
-	G4endl;
-      SetLowEnergyLimit(lowEnergyLimit);
-  }
- 
-  if (HighEnergyLimit() > highEnergyLimit)
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
+  }
+  
+  if (shellCrossSectionHandler)
   {
-      G4cout << "G4LivermorePhotoElectricModel: high energy limit decreased from " <<
-        HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" 
-	     << G4endl;
-      SetHighEnergyLimit(highEnergyLimit);
-  }
-
+    shellCrossSectionHandler->Clear();
+    delete shellCrossSectionHandler;
+  }
+  
   // Read data tables for all materials
   
@@ -105,7 +114,4 @@
   crossSectionHandler->LoadData(crossSectionFile);
 
-  meanFreePathTable = 0;
-  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
-
   shellCrossSectionHandler = new G4CrossSectionHandler();
   shellCrossSectionHandler->Clear();
@@ -113,26 +119,26 @@
   shellCrossSectionHandler->LoadShellData(shellCrossSectionFile);
   
-  // SI - Buggy default ?
+  // SI - Simple generator is buggy
   //generatorName = "geant4.6.2";
   //ElectronAngularGenerator = new G4PhotoElectricAngularGeneratorSimple("GEANTSimpleGenerator");              // default generator
-
-  //
-  
+  ElectronAngularGenerator = 
+    new G4PhotoElectricAngularGeneratorSauterGavrila("GEANTSauterGavrilaGenerator");        
+
+  //  
   if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for Livermore PhotoElectric model" << G4endl;
 
-  G4cout << "Livermore PhotoElectric model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  //  InitialiseElementSelectors(particle,cuts);
+
+  if (verboseLevel > 0) { 
+    G4cout << "Livermore PhotoElectric model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / eV << " eV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -147,5 +153,9 @@
 {
   if (verboseLevel > 3)
-    G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermorePhotoElectricModel" << G4endl;
+    G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermorePhotoElectricModel" 
+	   << G4endl;
+
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit)
+    return 0;
 
   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
@@ -155,9 +165,10 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-void G4LivermorePhotoElectricModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
-					      const G4MaterialCutsCouple* couple,
-					      const G4DynamicParticle* aDynamicGamma,
-					      G4double,
-					      G4double)
+void 
+G4LivermorePhotoElectricModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
+						 const G4MaterialCutsCouple* couple,
+						 const G4DynamicParticle* aDynamicGamma,
+						 G4double,
+						 G4double)
 {
 
@@ -172,17 +183,23 @@
   G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
   
+  // kill incident photon
+  fParticleChange->SetProposedKineticEnergy(0.);
+  fParticleChange->ProposeTrackStatus(fStopAndKill);   
+
+  // low-energy gamma is absorpted by this process
   if (photonEnergy <= lowEnergyLimit)
     {
-      fParticleChange->ProposeTrackStatus(fStopAndKill);
-      fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
-      // SI - IS THE FOLLOWING RETURN NECESSARY ?
-      return ;
+      return;
     }
  
-  G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection(); // Returns the normalized direction of the momentum
+  // Returns the normalized direction of the momentum
+  G4ThreeVector photonDirection = aDynamicGamma->GetMomentumDirection(); 
 
   // Select randomly one element in the current material
-  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
+  //  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
+  const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+  const G4Element* elm = SelectRandomAtom(couple->GetMaterial(),particle,photonEnergy);
+  G4int Z = (G4int)elm->GetZ();
 
   // Select the ionised shell in the current atom according to shell cross sections
@@ -195,11 +212,4 @@
   G4int shellId = shell->ShellId();
 
-  // Create lists of pointers to DynamicParticles (photons and electrons)
-  // (Is the electron vector necessary? To be checked)
-  std::vector<G4DynamicParticle*>* photonVector = 0;
-  std::vector<G4DynamicParticle*> electronVector;
-
-  G4double energyDeposit = 0.0;
-
   // Primary outcoming electron
   G4double eKineticEnergy = photonEnergy - bindingEnergy;
@@ -209,26 +219,17 @@
   if (eKineticEnergy > 0.)
     {
-      // SI - Removed safety
-      
-      // Generate the electron only if with large enough range w.r.t. cuts and safety
-      //G4double safety = aStep.GetPostStepPoint()->GetSafety();
-
-      //if (rangeTest->Escape(G4Electron::Electron(),couple,eKineticEnergy,safety))
-	{
-
-	  // Calculate direction of the photoelectron
-	  G4ThreeVector gammaPolarization = aDynamicGamma->GetPolarization();
-	  G4ThreeVector electronDirection = ElectronAngularGenerator->GetPhotoElectronDirection(photonDirection,eKineticEnergy,gammaPolarization,shellId);
-
-	  // The electron is created ...
-	  G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(),
-							       electronDirection,
-							       eKineticEnergy);
-	  electronVector.push_back(electron);
-	}
-      /*else
-	{
-	  energyDeposit += eKineticEnergy;
-	}*/
+      // Calculate direction of the photoelectron
+      G4ThreeVector gammaPolarization = aDynamicGamma->GetPolarization();
+      G4ThreeVector electronDirection = 
+	ElectronAngularGenerator->GetPhotoElectronDirection(photonDirection,
+							    eKineticEnergy,
+							    gammaPolarization,
+							    shellId);
+
+      // The electron is created ...
+      G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(),
+							   electronDirection,
+							   eKineticEnergy);
+      fvect->push_back(electron);
     }
   else
@@ -237,88 +238,54 @@
     }
 
-  G4int nElectrons = electronVector.size();
-  size_t nTotPhotons = 0;
-  G4int nPhotons=0;
-  const G4ProductionCutsTable* theCoupleTable=
-        G4ProductionCutsTable::GetProductionCutsTable();
-
-  size_t index = couple->GetIndex();
-  G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[index];
-  cutg = std::min(cutForLowEnergySecondaryPhotons,cutg);
-  
-  G4double cute = (*(theCoupleTable->GetEnergyCutsVector(1)))[index];
-  cute = std::min(cutForLowEnergySecondaryPhotons,cute);
-
-  G4DynamicParticle* aPhoton;
-
-  // Generation of fluorescence
-  // Data in EADL are available only for Z > 5
-  // Protection to avoid generating photons in the unphysical case of
-  // shell binding energy > photon energy
-  if (Z > 5  && (bindingEnergy > cutg || bindingEnergy > cute))
-    {
-      photonVector = deexcitationManager.GenerateParticles(Z,shellId);
-      nTotPhotons = photonVector->size();
-      for (size_t k=0; k<nTotPhotons; k++)
-	{
-	  aPhoton = (*photonVector)[k];
-	  if (aPhoton)
-	    {
-              G4double itsCut = cutg;
-              if(aPhoton->GetDefinition() == G4Electron::Electron()) itsCut = cute;
-	      G4double itsEnergy = aPhoton->GetKineticEnergy();
-
-	      if (itsEnergy > itsCut && itsEnergy <= bindingEnergy)
-		{
-		  nPhotons++;
-		  // Local energy deposit is given as the sum of the
-		  // energies of incident photons minus the energies
-		  // of the outcoming fluorescence photons
-		  bindingEnergy -= itsEnergy;
-
-		}
-	      else
-		{
-                  delete aPhoton;
-                  (*photonVector)[k] = 0;
-                }
-	    }
-	}
-    }
-
-  energyDeposit += bindingEnergy;
-
-  // Final state
-  
-  for (G4int l = 0; l<nElectrons; l++ )
-    {
-      aPhoton = electronVector[l];
-      if(aPhoton) {
-        fvect->push_back(aPhoton);
+  // deexcitation
+  if(DeexcitationFlag() && Z > 5) {
+
+    const G4ProductionCutsTable* theCoupleTable=
+      G4ProductionCutsTable::GetProductionCutsTable();
+
+    size_t index = couple->GetIndex();
+    G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[index];
+    G4double cute = (*(theCoupleTable->GetEnergyCutsVector(1)))[index];
+
+    // Generation of fluorescence
+    // Data in EADL are available only for Z > 5
+    // Protection to avoid generating photons in the unphysical case of
+    // shell binding energy > photon energy
+    if (bindingEnergy > cutg || bindingEnergy > cute)
+      {
+	G4DynamicParticle* aPhoton;
+	deexcitationManager.SetCutForSecondaryPhotons(cutg);
+	deexcitationManager.SetCutForAugerElectrons(cute);
+ 
+	std::vector<G4DynamicParticle*>* photonVector = 
+	  deexcitationManager.GenerateParticles(Z,shellId);
+	size_t nTotPhotons = photonVector->size();
+	for (size_t k=0; k<nTotPhotons; k++)
+	  {
+	    aPhoton = (*photonVector)[k];
+	    if (aPhoton)
+	      {
+		G4double itsEnergy = aPhoton->GetKineticEnergy();
+		if (itsEnergy <= bindingEnergy)
+		  {
+		    // Local energy deposit is given as the sum of the
+		    // energies of incident photons minus the energies
+		    // of the outcoming fluorescence photons
+		    bindingEnergy -= itsEnergy;
+		    fvect->push_back(aPhoton);
+		  }
+		else
+		  {
+		    // abnormal case of energy non-conservation
+		    delete aPhoton;
+		    (*photonVector)[k] = 0;
+		  }
+	      }
+	  }
+	delete photonVector;
       }
-    }
-  for ( size_t ll = 0; ll < nTotPhotons; ll++)
-    {
-      aPhoton = (*photonVector)[ll];
-      if(aPhoton) {
-        fvect->push_back(aPhoton);
-      }
-    }
-
-  delete photonVector;
-
-  if (energyDeposit < 0)
-    {
-      G4cout << "WARNING - "
-	     << "G4LowEnergyPhotoElectric::PostStepDoIt - Negative energy deposit"
-	     << G4endl;
-      energyDeposit = 0;
-    }
-
-  // kill incident photon
-  fParticleChange->ProposeMomentumDirection( 0., 0., 0. );
-  fParticleChange->SetProposedKineticEnergy(0.);
-  fParticleChange->ProposeTrackStatus(fStopAndKill);   
-  fParticleChange->ProposeLocalEnergyDeposit(energyDeposit);
+  }
+  // excitation energy left
+  fParticleChange->ProposeLocalEnergyDeposit(bindingEnergy);
 }
 
@@ -348,7 +315,8 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-void G4LivermorePhotoElectricModel::SetAngularGenerator(G4VPhotoElectricAngularDistribution* distribution)
-{
-  ElectronAngularGenerator = distribution;
+void 
+G4LivermorePhotoElectricModel::SetAngularGenerator(G4VPhotoElectricAngularDistribution* dist)
+{
+  ElectronAngularGenerator = dist;
   ElectronAngularGenerator->PrintGeneratorInformation();
 }
@@ -361,5 +329,6 @@
     {
       delete ElectronAngularGenerator;
-      ElectronAngularGenerator = new G4PhotoElectricAngularGeneratorSimple("GEANT4LowEnergySimpleGenerator");
+      ElectronAngularGenerator = 
+	new G4PhotoElectricAngularGeneratorSimple("GEANT4LowEnergySimpleGenerator");
       generatorName = name;
     }
@@ -367,5 +336,6 @@
     {
       delete ElectronAngularGenerator;
-      ElectronAngularGenerator = new G4PhotoElectricAngularGeneratorSauterGavrila("GEANT4SauterGavrilaGenerator");
+      ElectronAngularGenerator = 
+	new G4PhotoElectricAngularGeneratorSauterGavrila("GEANT4SauterGavrilaGenerator");
       generatorName = name;
     }
@@ -373,5 +343,6 @@
     {
       delete ElectronAngularGenerator;
-      ElectronAngularGenerator = new G4PhotoElectricAngularGeneratorPolarized("GEANT4LowEnergyPolarizedGenerator");
+      ElectronAngularGenerator = 
+	new G4PhotoElectricAngularGeneratorPolarized("GEANT4LowEnergyPolarizedGenerator");
       generatorName = name;
     }
@@ -384,26 +355,2 @@
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermorePhotoElectricModel::GetMeanFreePath(const G4Track& track,
-						   G4double, // previousStepSize
-					       G4ForceCondition*)
-{
-  const G4DynamicParticle* photon = track.GetDynamicParticle();
-  G4double energy = photon->GetKineticEnergy();
-  G4Material* material = track.GetMaterial();
-  //  size_t materialIndex = material->GetIndex();
-
-  G4double meanFreePath = DBL_MAX;
-
-  //  if (energy > highEnergyLimit) 
-  //    meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
-  //  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
-  //  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
-
-  G4double cross = shellCrossSectionHandler->ValueForMaterial(material,energy);
-  if(cross > 0.0) meanFreePath = 1.0/cross;
-
-  return meanFreePath;
-}
-
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedComptonModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedComptonModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedComptonModel.cc	(revision 1055)
@@ -24,7 +24,17 @@
 // ********************************************************************
 //
-// $Id: G4LivermorePolarizedComptonModel.cc,v 1.1 2008/10/30 14:16:35 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LivermorePolarizedComptonModel.cc,v 1.6 2009/05/03 08:29:55 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
+// History:
+// --------
+// 02 May 2009   S Incerti as V. Ivanchenko proposed in G4LivermoreComptonModel.cc
+//
+// Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - remove GetMeanFreePath method and table
+//                  - added protection against numerical problem in energy sampling 
+//                  - use G4ElementSelector
 
 #include "G4LivermorePolarizedComptonModel.hh"
@@ -38,9 +48,9 @@
 G4LivermorePolarizedComptonModel::G4LivermorePolarizedComptonModel(const G4ParticleDefinition*,
                                              const G4String& nam)
-:G4VEmModel(nam),isInitialised(false)
-{
-  lowEnergyLimit = 250 * eV; // SI - Could be 10 eV ?
+:G4VEmModel(nam),isInitialised(false),meanFreePathTable(0),scatterFunctionData(0),crossSectionHandler(0)
+{
+  lowEnergyLimit = 250 * eV; 
   highEnergyLimit = 100 * GeV;
-  SetLowEnergyLimit(lowEnergyLimit);
+  //SetLowEnergyLimit(lowEnergyLimit);
   SetHighEnergyLimit(highEnergyLimit);
   
@@ -53,10 +63,11 @@
   // 4 = entering in methods
 
+  if( verboseLevel>0 ) { 
   G4cout << "Livermore Polarized Compton is constructed " << G4endl
          << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
+         << lowEnergyLimit / eV << " eV - "
          << highEnergyLimit / GeV << " GeV"
          << G4endl;
-
+  }
 }
 
@@ -65,7 +76,7 @@
 G4LivermorePolarizedComptonModel::~G4LivermorePolarizedComptonModel()
 {  
-  delete meanFreePathTable;
-  delete crossSectionHandler;
-  delete scatterFunctionData;
+  if (meanFreePathTable)   delete meanFreePathTable;
+  if (crossSectionHandler) delete crossSectionHandler;
+  if (scatterFunctionData) delete scatterFunctionData;
 }
 
@@ -78,21 +89,9 @@
     G4cout << "Calling G4LivermorePolarizedComptonModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
-  
-  if (LowEnergyLimit() < lowEnergyLimit)
-    {
-      G4cout << "G4LivermorePolarizedComptonModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
-      SetLowEnergyLimit(lowEnergyLimit);
-    }
-
-  if (HighEnergyLimit() > highEnergyLimit)
-    {
-      G4cout << "G4LivermorePolarizedComptonModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
-      SetHighEnergyLimit(highEnergyLimit);
-    }
+  if (crossSectionHandler)
+  {
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
+  }
 
   // Reading of data files - all materials are read
@@ -116,23 +115,21 @@
   shellData.LoadData(file);
 
-  //
   if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for Livermore Polarized Compton model" << G4endl;
 
-  G4cout << "Livermore Polarized Compton model is initialized " << G4endl
+  InitialiseElementSelectors(particle,cuts);
+
+  if(  verboseLevel>0 ) { 
+    G4cout << "Livermore Polarized Compton model is initialized " << G4endl
          << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
+         << LowEnergyLimit() / eV << " eV - "
          << HighEnergyLimit() / GeV << " GeV"
          << G4endl;
-
+  }
+  
   //
     
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-    
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -149,4 +146,6 @@
     G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermorePolarizedComptonModel" << G4endl;
 
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit) return 0.0;
+
   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
   return cs;
@@ -203,5 +202,4 @@
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy0);
-      // SI - IS THE FOLLOWING RETURN NECESSARY ?
       return;
     }
@@ -210,6 +208,8 @@
 
   // Select randomly one element in the current material
-
-  G4int Z = crossSectionHandler->SelectRandomAtom(couple,gammaEnergy0);
+  //G4int Z = crossSectionHandler->SelectRandomAtom(couple,gammaEnergy0);
+  const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+  const G4Element* elm = SelectRandomAtom(couple,particle,gammaEnergy0);
+  G4int Z = (G4int)elm->GetZ();
 
   // Sample the energy and the polarization of the scattered photon
@@ -429,4 +429,5 @@
   else
     {
+      gammaEnergy1 = 0.;
       fParticleChange->SetProposedKineticEnergy(0.) ;
       fParticleChange->ProposeTrackStatus(fStopAndKill);
@@ -439,4 +440,11 @@
   G4double ElecKineEnergy = gammaEnergy0 - gammaEnergy1 -bindingE;
 
+  // SI -protection against negative final energy: no e- is created
+  // like in G4LivermoreComptonModel.cc
+  if(ElecKineEnergy < 0.0) {
+    fParticleChange->ProposeLocalEnergyDeposit(gammaEnergy0 - gammaEnergy1);
+    return;
+  }
+ 
   // SI - Removed range test
   
@@ -664,20 +672,3 @@
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermorePolarizedComptonModel::GetMeanFreePath(const G4Track& track,
-						      G4double,
-						      G4ForceCondition*)
-{
-  const G4DynamicParticle* photon = track.GetDynamicParticle();
-  G4double energy = photon->GetKineticEnergy();
-  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
-  size_t materialIndex = couple->GetIndex();
-  G4double meanFreePath;
-  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
-  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
-  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
-  return meanFreePath;
-}
-
-
+
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedRayleighModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedRayleighModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermorePolarizedRayleighModel.cc	(revision 1055)
@@ -24,7 +24,17 @@
 // ********************************************************************
 //
-// $Id: G4LivermorePolarizedRayleighModel.cc,v 1.1 2008/10/30 14:16:35 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LivermorePolarizedRayleighModel.cc,v 1.5 2009/05/02 15:20:53 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
+// History:
+// --------
+// 02 May 2009   S Incerti as V. Ivanchenko proposed in G4LivermoreRayleighModel.cc
+//
+// Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - remove GetMeanFreePath method and table
+//                  - remove initialisation of element selector 
+//                  - use G4ElementSelector
 
 #include "G4LivermorePolarizedRayleighModel.hh"
@@ -38,10 +48,10 @@
 G4LivermorePolarizedRayleighModel::G4LivermorePolarizedRayleighModel(const G4ParticleDefinition*,
                                              const G4String& nam)
-:G4VEmModel(nam),isInitialised(false)
-{
-  lowEnergyLimit = 250 * eV; // SI - Could be 10 eV ?
+:G4VEmModel(nam),isInitialised(false),crossSectionHandler(0),formFactorData(0)
+{
+  lowEnergyLimit = 250 * eV; 
   highEnergyLimit = 100 * GeV;
   
-  SetLowEnergyLimit(lowEnergyLimit);
+  //SetLowEnergyLimit(lowEnergyLimit);
   SetHighEnergyLimit(highEnergyLimit);
   //
@@ -54,9 +64,11 @@
   // 4 = entering in methods
 
-  G4cout << "Livermore Polarized Rayleigh is constructed " << G4endl
+  if(verboseLevel > 0) {
+    G4cout << "Livermore Polarized Rayleigh is constructed " << G4endl
          << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
+         << lowEnergyLimit / eV << " eV - "
          << highEnergyLimit / GeV << " GeV"
          << G4endl;
+  }
 }
 
@@ -65,6 +77,6 @@
 G4LivermorePolarizedRayleighModel::~G4LivermorePolarizedRayleighModel()
 {  
-  delete crossSectionHandler;
-  delete formFactorData;
+  if (crossSectionHandler) delete crossSectionHandler;
+  if (formFactorData) delete formFactorData;
 }
 
@@ -84,21 +96,10 @@
     G4cout << "Calling G4LivermorePolarizedRayleighModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
-  
-  if (LowEnergyLimit() < lowEnergyLimit)
-    {
-      G4cout << "G4LivermorePolarizedRayleighModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
-      SetLowEnergyLimit(lowEnergyLimit);
-    }
-  if (HighEnergyLimit() > highEnergyLimit)
-    {
-      G4cout << "G4LivermorePolarizedRayleighModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
-      SetHighEnergyLimit(highEnergyLimit);
-    }
-    
+  if (crossSectionHandler)
+  {
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
+  }
+  
   // Read data files for all materials
 
@@ -113,21 +114,22 @@
   formFactorData->LoadData(formFactorFile);
 
+  InitialiseElementSelectors(particle,cuts);
+
   //
   if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for Livermore Polarized Rayleigh model" << G4endl;
 
-  G4cout << "Livermore Polarized Rayleigh model is initialized " << G4endl
+  InitialiseElementSelectors(particle,cuts);
+
+  if (verboseLevel > 0) { 
+    G4cout << "Livermore Polarized Rayleigh model is initialized " << G4endl
          << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
+         << LowEnergyLimit() / eV << " eV - "
          << HighEnergyLimit() / GeV << " GeV"
          << G4endl;
+	 }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-    
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -144,4 +146,6 @@
     G4cout << "Calling CrossSectionPerAtom() of G4LivermorePolarizedRayleighModel" << G4endl;
 
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit) return 0.0;
+
   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
   return cs;
@@ -166,5 +170,4 @@
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
-      // SI - IS THE FOLLOWING RETURN NECESSARY ?
       return ;
   }
@@ -173,5 +176,8 @@
 
   // Select randomly one element in the current material
-  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  // G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+  const G4Element* elm = SelectRandomAtom(couple,particle,photonEnergy0);
+  G4int Z = (G4int)elm->GetZ();
 
   G4double outcomingPhotonCosTheta = GenerateCosTheta(photonEnergy0, Z);
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreRayleighModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreRayleighModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LivermoreRayleighModel.cc	(revision 1055)
@@ -24,7 +24,18 @@
 // ********************************************************************
 //
-// $Id: G4LivermoreRayleighModel.cc,v 1.1 2008/10/30 14:16:35 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
-//
+// $Id: G4LivermoreRayleighModel.cc,v 1.6 2009/04/18 18:29:34 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
+//
+// Author: Sebastien Inserti
+//         30 October 2008
+//
+// History:
+// --------
+// 18 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - remove GetMeanFreePath method and table
+//                  - remove initialisation of element selector 
+//                  - use G4ElementSelector
 
 #include "G4LivermoreRayleighModel.hh"
@@ -37,11 +48,12 @@
 
 G4LivermoreRayleighModel::G4LivermoreRayleighModel(const G4ParticleDefinition*,
-                                             const G4String& nam)
-:G4VEmModel(nam),isInitialised(false)
-{
-  lowEnergyLimit = 250 * eV; // SI - Could be 10 eV ?
+						   const G4String& nam)
+  :G4VEmModel(nam),isInitialised(false),meanFreePathTable(0),
+   formFactorData(0),crossSectionHandler(0)
+{
+  lowEnergyLimit = 250 * eV; 
   highEnergyLimit = 100 * GeV;
   
-  SetLowEnergyLimit(lowEnergyLimit);
+  //  SetLowEnergyLimit(lowEnergyLimit);
   SetHighEnergyLimit(highEnergyLimit);
   //
@@ -54,9 +66,11 @@
   // 4 = entering in methods
 
-  G4cout << "Livermore Rayleigh is constructed " << G4endl
-         << "Energy range: "
-         << lowEnergyLimit / keV << " keV - "
-         << highEnergyLimit / GeV << " GeV"
-         << G4endl;
+  if(verboseLevel > 0) {
+    G4cout << "Livermore Rayleigh is constructed " << G4endl
+	   << "Energy range: "
+	   << lowEnergyLimit / eV << " eV - "
+	   << highEnergyLimit / GeV << " GeV"
+	   << G4endl;
+  }
 }
 
@@ -65,7 +79,6 @@
 G4LivermoreRayleighModel::~G4LivermoreRayleighModel()
 {  
-  delete meanFreePathTable;
-  delete crossSectionHandler;
-  delete formFactorData;
+  if (crossSectionHandler) delete crossSectionHandler;
+  if (formFactorData) delete formFactorData;
 }
 
@@ -73,27 +86,15 @@
 
 void G4LivermoreRayleighModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+					  const G4DataVector& cuts)
 {
   if (verboseLevel > 3)
     G4cout << "Calling G4LivermoreRayleighModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
-  // Energy limits
+  if (crossSectionHandler)
+  {
+    crossSectionHandler->Clear();
+    delete crossSectionHandler;
+  }
   
-  if (LowEnergyLimit() < lowEnergyLimit)
-  {
-    G4cout << "G4LivermoreRayleighModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
-    SetLowEnergyLimit(lowEnergyLimit);
-  }
-
-  if (HighEnergyLimit() > highEnergyLimit)
-  {
-    G4cout << "G4LivermoreRayleighModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
-    SetHighEnergyLimit(highEnergyLimit);
-  }
-
   // Data are read for all materials
   
@@ -103,7 +104,4 @@
   crossSectionHandler->LoadData(crossSectionFile);
 
-  meanFreePathTable = 0;
-  meanFreePathTable = crossSectionHandler->BuildMeanFreePathForMaterials();
-
   G4VDataSetAlgorithm* ffInterpolation = new G4LogLogInterpolation;
   G4String formFactorFile = "rayl/re-ff-";
@@ -111,22 +109,20 @@
   formFactorData->LoadData(formFactorFile);
 
-  //
-  
+  InitialiseElementSelectors(particle,cuts);
+
+  //  
   if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for Livermore Rayleigh model" << G4endl;
 
-  G4cout << "Livermore Rayleigh model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if (verboseLevel > 0) { 
+    G4cout << "Livermore Rayleigh model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / eV << " eV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
-
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 
@@ -144,4 +140,7 @@
     G4cout << "Calling CrossSectionPerAtom() of G4LivermoreRayleighModel" << G4endl;
 
+  if (GammaEnergy < lowEnergyLimit || GammaEnergy > highEnergyLimit)
+    return 0.0;
+
   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
   return cs;
@@ -160,18 +159,21 @@
 
   G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
-  
+
+  // absorption of low-energy gamma  
   if (photonEnergy0 <= lowEnergyLimit)
-  {
+    {
       fParticleChange->ProposeTrackStatus(fStopAndKill);
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
-      // SI - IS THE FOLLOWING RETURN NECESSARY ?
       return ;
-  }
+    }
 
   G4ParticleMomentum photonDirection0 = aDynamicGamma->GetMomentumDirection();
 
   // Select randomly one element in the current material
-  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  //  G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy0);
+  const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
+  const G4Element* elm = SelectRandomAtom(couple,particle,photonEnergy0);
+  G4int Z = (G4int)elm->GetZ();
 
   // Sample the angle of the scattered photon
@@ -188,15 +190,19 @@
     {
       do
-      {
-      cosTheta = 2. * G4UniformRand() - 1.;
-      fcostheta = ( 1. + cosTheta*cosTheta)/2.;
-      } while (fcostheta < G4UniformRand());
+	{
+	  cosTheta = 2. * G4UniformRand() - 1.;
+	  fcostheta = ( 1. + cosTheta*cosTheta)/2.;
+	} while (fcostheta < G4UniformRand());
 
       G4double sinThetaHalf = std::sqrt((1. - cosTheta) / 2.);
       x = sinThetaHalf / (wlPhoton/cm);
       if (x > 1.e+005)
-         dataFormFactor = formFactorData->FindValue(x,Z-1);
+	{
+	  dataFormFactor = formFactorData->FindValue(x,Z-1);
+	}
       else
-         dataFormFactor = formFactorData->FindValue(0.,Z-1);
+	{
+	  dataFormFactor = formFactorData->FindValue(0.,Z-1);
+	}
       randomFormFactor = G4UniformRand() * Z * Z;
       sinTheta = std::sqrt(1. - cosTheta*cosTheta);
@@ -219,20 +225,3 @@
 }
 
-//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
-
-G4double G4LivermoreRayleighModel::GetMeanFreePath(const G4Track& track, 
-					      G4double, // previousStepSize
-					      G4ForceCondition*)
-{
-  const G4DynamicParticle* photon = track.GetDynamicParticle();
-  G4double energy = photon->GetKineticEnergy();
-  const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
-  size_t materialIndex = couple->GetIndex();
-
-  G4double meanFreePath;
-  if (energy > highEnergyLimit) meanFreePath = meanFreePathTable->FindValue(highEnergyLimit,materialIndex);
-  else if (energy < lowEnergyLimit) meanFreePath = DBL_MAX;
-  else meanFreePath = meanFreePathTable->FindValue(energy,materialIndex);
-  return meanFreePath;
-}
-
+
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyBremsstrahlung.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyBremsstrahlung.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyBremsstrahlung.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4LowEnergyBremsstrahlung.cc,v 1.71 2006/06/29 19:40:13 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyBremsstrahlung.cc,v 1.72 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // --------------------------------------------------------------
@@ -93,4 +93,15 @@
 //  angularDistribution->PrintGeneratorInformation();
   TsaiAngularDistribution = new G4ModifiedTsai("TsaiGenerator");
+
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyBremsstrahlung is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyCompton.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyCompton.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyCompton.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4LowEnergyCompton.cc,v 1.47 2008/12/18 13:01:28 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyCompton.cc,v 1.48 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: A. Forti
@@ -104,4 +104,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyCompton is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyGammaConversion.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyGammaConversion.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyGammaConversion.cc	(revision 1055)
@@ -26,6 +26,6 @@
 // --------------------------------------------------------------------
 ///
-// $Id: G4LowEnergyGammaConversion.cc,v 1.36 2006/06/29 19:40:17 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyGammaConversion.cc,v 1.37 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // 
@@ -104,4 +104,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyGammaConversion is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
  
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyIonisation.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyIonisation.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyIonisation.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4LowEnergyIonisation.cc,v 1.103 2008/05/02 19:23:38 pia Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyIonisation.cc,v 1.104 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // --------------------------------------------------------------
@@ -131,4 +131,14 @@
   cutForElectrons = 250.0*eV;
   verboseLevel = 0;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyIonisation is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPhotoElectric.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPhotoElectric.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPhotoElectric.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4LowEnergyPhotoElectric.cc,v 1.56 2006/06/29 19:40:23 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyPhotoElectric.cc,v 1.57 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: A. Forti
@@ -122,4 +122,14 @@
 	     << G4endl;
     }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyPhotoElectric is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedCompton.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedCompton.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedCompton.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4LowEnergyPolarizedCompton.cc,v 1.25 2008/05/02 19:23:38 pia Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyPolarizedCompton.cc,v 1.26 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // ------------------------------------------------------------
@@ -122,4 +122,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyPolarizedCompton is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedRayleigh.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedRayleigh.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyPolarizedRayleigh.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4LowEnergyPolarizedRayleigh.cc,v 1.7 2006/06/29 19:40:27 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyPolarizedRayleigh.cc,v 1.8 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // --------------------------------------------------------------
@@ -63,4 +63,13 @@
  
  
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyPolarizedRayleigh is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyRayleigh.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyRayleigh.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4LowEnergyRayleigh.cc	(revision 1055)
@@ -26,6 +26,6 @@
 // --------------------------------------------------------------------
 //
-// $Id: G4LowEnergyRayleigh.cc,v 1.37 2006/06/29 19:40:29 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4LowEnergyRayleigh.cc,v 1.38 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: A. Forti
@@ -99,4 +99,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4LowEnergyRayleigh is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilation.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilation.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilation.cc	(revision 1055)
@@ -64,4 +64,14 @@
 	     << G4endl;
     }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeAnnihilation is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilationModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilationModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeAnnihilationModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeAnnihilationModel.cc,v 1.2 2008/12/04 14:09:36 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeAnnihilationModel.cc,v 1.3 2009/04/17 10:29:20 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -32,5 +32,8 @@
 // --------
 // 29 Oct 2008   L Pandola    Migration from process to model 
-//
+// 15 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - do not use G4ElementSelector
 
 #include "G4PenelopeAnnihilationModel.hh"
@@ -48,7 +51,7 @@
   :G4VEmModel(nam),isInitialised(false)
 {
-  fIntrinsicLowEnergyLimit = 0.0*eV;
+  fIntrinsicLowEnergyLimit = 0.0;
   fIntrinsicHighEnergyLimit = 100.0*GeV;
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
  
@@ -73,36 +76,20 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-void G4PenelopeAnnihilationModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+void G4PenelopeAnnihilationModel::Initialise(const G4ParticleDefinition*,
+					     const G4DataVector&)
 {
   if (verboseLevel > 3)
     G4cout << "Calling G4PenelopeAnnihilationModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopeAnnihilationModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopeAnnihilationModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
-
-  G4cout << "Penelope Annihilation model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if(verboseLevel > 0) {
+    G4cout << "Penelope Annihilation model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / keV << " keV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true; 
 }
@@ -155,6 +142,10 @@
 
   G4double kineticEnergy = aDynamicPositron->GetKineticEnergy();
-  
-  if (kineticEnergy == 0)
+
+  // kill primary
+  fParticleChange->SetProposedKineticEnergy(0.);
+  fParticleChange->ProposeTrackStatus(fStopAndKill);
+  
+  if (kineticEnergy == 0.0)
     {
       //Old AtRestDoIt
@@ -170,6 +161,4 @@
       fvect->push_back(firstGamma);
       fvect->push_back(secondGamma);
-      fParticleChange->SetProposedKineticEnergy(0.);
-      fParticleChange->ProposeTrackStatus(fStopAndKill);
       return;
     }
@@ -227,6 +216,4 @@
 							   photon2Energy);
   fvect->push_back(aParticle2);
-  fParticleChange->SetProposedKineticEnergy(0.);
-  fParticleChange->ProposeTrackStatus(fStopAndKill);
 
   if (verboseLevel > 1)
@@ -244,6 +231,5 @@
     }
   if (verboseLevel > 0)
-    {
-      
+    {      
       G4double energyDiff = std::fabs(totalAvailableEnergy-photon1Energy-photon2Energy);
       if (energyDiff > 0.05*keV)
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlung.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlung.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlung.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeBremsstrahlung.cc,v 1.18 2006/06/29 19:40:35 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeBremsstrahlung.cc,v 1.19 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // --------------------------------------------------------------
@@ -71,4 +71,15 @@
   cutForPhotons = 0.;
   verboseLevel = 0;
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeBremsstrahlung is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
+
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlungContinuous.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlungContinuous.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeBremsstrahlungContinuous.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeBremsstrahlungContinuous.cc,v 1.10 2008/12/09 15:08:13 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeBremsstrahlungContinuous.cc,v 1.11 2008/12/15 09:23:06 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // --------------------------------------------------------------
@@ -97,6 +97,5 @@
  G4String dirFile = pathString + "/penelope/" + name;
  std::ifstream file(dirFile);
- std::filebuf* lsdp = file.rdbuf();
- if (!(lsdp->is_open()))
+ if (!file.is_open())
      {
       G4String excep = "G4PenelopeBremsstrahlungContinuous - data file " + name + " not found!";
@@ -112,14 +111,15 @@
      for (size_t j=0;j<NumberofKPoints;j++){
        file >> a1;
-       ReducedCS[i][j]=a1*cm2; //coversion present in Penelope source
+       ReducedCS[i][j]=a1*cm2; 
      }
      //3) read the total cross section, in cm2
      file >> a1;
-     TotalCS[i]=a1*cm2; //conversion present in Penelope source
+     TotalCS[i]=a1*cm2; 
      // Check closing item
      file >> a1;
      if (a1 != ((G4double) -1))
        {
-	 G4String excep = "G4PenelopeBremsstrahlungContinuous - Check the bremms data file " + name;
+	 G4String excep = "G4PenelopeBremsstrahlungContinuous - Check the bremms data file " 
+	   + name;
 	 G4Exception(excep);
        }
@@ -229,5 +229,6 @@
   
   //Global x-section factor
-  G4double Fact=Zmat*Zmat*(energy+electron_mass_c2)*(energy+electron_mass_c2)/(energy*(energy+2.0*electron_mass_c2));
+  G4double Fact=Zmat*Zmat*(energy+electron_mass_c2)*(energy+electron_mass_c2)/
+    (energy*(energy+2.0*electron_mass_c2));
   Fact *= PositronCorrection(energy);
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeCompton.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeCompton.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeCompton.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeCompton.cc,v 1.33 2008/06/03 15:44:25 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeCompton.cc,v 1.34 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -120,4 +120,16 @@
 	     << G4endl;
     }
+
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeCompton is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
+
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeComptonModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeComptonModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeComptonModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeComptonModel.cc,v 1.2 2008/12/04 14:11:21 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeComptonModel.cc,v 1.4 2009/04/18 18:29:34 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -40,4 +40,8 @@
 //			      that are read from the external files.
 // 24 Nov 2008   L Pandola    Find a cleaner way to delete vectors.
+// 16 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - do not apply production threshold on level of the model
 //
 
@@ -69,5 +73,5 @@
   fIntrinsicLowEnergyLimit = 100.0*eV;
   fIntrinsicHighEnergyLimit = 100.0*GeV;
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
   //
@@ -118,5 +122,5 @@
 
 void G4PenelopeComptonModel::Initialise(const G4ParticleDefinition* particle,
-                                       const G4DataVector& cuts)
+					const G4DataVector& cuts)
 {
   if (verboseLevel > 3)
@@ -124,29 +128,15 @@
 
   InitialiseElementSelectors(particle,cuts);
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopeComptonModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopeComptonModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
-
-  G4cout << "Penelope Compton model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+
+  if (verboseLevel > 0) {
+    G4cout << "Penelope Compton model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / keV << " keV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true; 
 }
@@ -249,11 +239,11 @@
   G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
 
-  if (photonEnergy0 <= LowEnergyLimit())
-  {
+  if (photonEnergy0 <= fIntrinsicLowEnergyLimit)
+    {
       fParticleChange->ProposeTrackStatus(fStopAndKill);
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
       return ;
-  }
+    }
 
   G4ParticleMomentum photonDirection0 = aDynamicGamma->GetMomentumDirection();
@@ -263,5 +253,6 @@
     G4cout << "Going to select element in " << couple->GetMaterial()->GetName() << G4endl;
   // atom can be selected effitiantly if element selectors are initialised 
-  const G4Element* anElement = SelectRandomAtom(couple,G4Gamma::GammaDefinition(),photonEnergy0);
+  const G4Element* anElement = 
+    SelectRandomAtom(couple,G4Gamma::GammaDefinition(),photonEnergy0);
   G4int Z = (G4int) anElement->GetZ();
   if (verboseLevel > 2)
@@ -320,5 +311,5 @@
       }while((epsilon*photonEnergy0-photonEnergy0+ionEnergy) >0);
     }
-  else //photonEnergy0<5 MeV
+  else //photonEnergy0 < 5 MeV
     {
       //Incoherent scattering function for theta=PI
@@ -472,9 +463,9 @@
   }
   
-
   //Create scattered electron
   G4double diffEnergy = photonEnergy0*(1-epsilon);
   ionEnergy = (*(ionizationEnergy->find(Z)->second))[iosc];
-  G4double Q2 = photonEnergy0*photonEnergy0+photonEnergy1*(photonEnergy1-2.0*photonEnergy0*cosTheta);
+  G4double Q2 = 
+    photonEnergy0*photonEnergy0+photonEnergy1*(photonEnergy1-2.0*photonEnergy0*cosTheta);
   G4double cosThetaE; //scattering angle for the electron
   if (Q2 > 1.0e-12)
@@ -496,8 +487,10 @@
   G4int shellId = shell->ShellId();
   G4double ionEnergyInPenelopeDatabase = ionEnergy;
-  ionEnergy = std::max(bindingEnergy,ionEnergyInPenelopeDatabase); //protection against energy non-conservation 
-
-  G4double eKineticEnergy = diffEnergy - ionEnergy; //subtract the excitation energy. If not emitted by fluorescence,
+  //protection against energy non-conservation
+  ionEnergy = std::max(bindingEnergy,ionEnergyInPenelopeDatabase);  
+
+  //subtract the excitation energy. If not emitted by fluorescence
   //the ionization energy is deposited as local energy deposition
+  G4double eKineticEnergy = diffEnergy - ionEnergy; 
   G4double localEnergyDeposit = ionEnergy; 
   G4double energyInFluorescence = 0.; //testing purposes only
@@ -505,5 +498,6 @@
   if (eKineticEnergy < 0) 
     {
-      //It means that there was some problem/mismatch between the two databases. Try to make it work
+      //It means that there was some problem/mismatch between the two databases. 
+      //Try to make it work
       //In this case available Energy (diffEnergy) < ionEnergy
       //Full residual energy is deposited locally
@@ -511,54 +505,61 @@
       eKineticEnergy = 0.0;
     }
-  G4double cutForLowEnergySecondaryPhotons = 250.0*eV;
-
-  //Get the cut for electrons
-  const G4ProductionCutsTable* theCoupleTable=
-    G4ProductionCutsTable::GetProductionCutsTable();
-  size_t indx = couple->GetIndex();
-  G4double cutE = (*(theCoupleTable->GetEnergyCutsVector(1)))[indx];
-  cutE = std::max(cutForLowEnergySecondaryPhotons,cutE);
      
   //the local energy deposit is what remains: part of this may be spent for fluorescence.
-  if (fUseAtomicDeexcitation)
-    { 
-      G4int nPhotons=0;
-
-      G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[indx];
-      cutg = std::max(cutForLowEnergySecondaryPhotons,cutg);
-
-      G4DynamicParticle* aPhoton;
-      G4AtomicDeexcitation deexcitationManager;
+  if(DeexcitationFlag() && Z > 5) {
+
+    const G4ProductionCutsTable* theCoupleTable=
+      G4ProductionCutsTable::GetProductionCutsTable();
+
+    size_t index = couple->GetIndex();
+    G4double cutg = (*(theCoupleTable->GetEnergyCutsVector(0)))[index];
+    G4double cute = (*(theCoupleTable->GetEnergyCutsVector(1)))[index];
+
+    // Generation of fluorescence
+    // Data in EADL are available only for Z > 5
+    // Protection to avoid generating photons in the unphysical case of
+    // shell binding energy > photon energy
+    if (localEnergyDeposit > cutg || localEnergyDeposit > cute)
+      { 
+	G4DynamicParticle* aPhoton;
+	G4AtomicDeexcitation deexcitationManager;
+	deexcitationManager.SetCutForSecondaryPhotons(cutg);
+	deexcitationManager.SetCutForAugerElectrons(cute);
+	deexcitationManager.ActivateAugerElectronProduction(false);
       
-      if (Z>5 && (localEnergyDeposit > cutg || localEnergyDeposit > cutE))
-	{
-	  photonVector = deexcitationManager.GenerateParticles(Z,shellId);
-	  for (size_t k=0;k<photonVector->size();k++){
-	    aPhoton = (*photonVector)[k];
-	    if (aPhoton)
+	photonVector = deexcitationManager.GenerateParticles(Z,shellId);
+	if(photonVector) 
+	  {
+	    size_t nPhotons = photonVector->size();
+	    for (size_t k=0; k<nPhotons; k++)
 	      {
-		G4double itsCut = cutg;
-		if (aPhoton->GetDefinition() == G4Electron::Electron()) itsCut = cutE;
-		G4double itsEnergy = aPhoton->GetKineticEnergy();
-		if (itsEnergy > itsCut && itsEnergy <= localEnergyDeposit)
+		aPhoton = (*photonVector)[k];
+		if (aPhoton)
 		  {
-		    nPhotons++;
-		    localEnergyDeposit -= itsEnergy;
-		    energyInFluorescence += itsEnergy;
-		  }
-		else
-		  {
-		    delete aPhoton;
-		    (*photonVector)[k]=0;
+		    G4double itsEnergy = aPhoton->GetKineticEnergy();
+		    if (itsEnergy <= localEnergyDeposit)
+		      {
+			localEnergyDeposit -= itsEnergy;
+			if (aPhoton->GetDefinition() == G4Gamma::Gamma()) 
+			  energyInFluorescence += itsEnergy;;
+			fvect->push_back(aPhoton);		    
+		      }
+		    else
+		      {
+			delete aPhoton;
+			(*photonVector)[k]=0;
+		      }
 		  }
 	      }
-	  }
-	}
-    }
+	    delete photonVector;
+	  }
+      }
+  }
 
   //Produce explicitely the electron only if its proposed kinetic energy is 
   //above the cut, otherwise add local energy deposition
   G4DynamicParticle* electron = 0;
-  if (eKineticEnergy > cutE)
+  //  if (eKineticEnergy > cutE) // VI: may be consistency problem if cut is applied here
+  if (eKineticEnergy > 0.0)
     {
       G4double xEl = sinThetaE * std::cos(phi+pi); 
@@ -576,18 +577,4 @@
     }
 
-  //This block below is executed only if there is at least one secondary photon produced by
-  //AtomicDeexcitation
-  if (photonVector)
-    {
-      for (size_t ll=0;ll<photonVector->size();ll++)
-	{
-	  if ((*photonVector)[ll]) 
-	    {
-	      G4DynamicParticle* aFluorescencePhoton = (*photonVector)[ll];
-	      fvect->push_back(aFluorescencePhoton);
-	    }
-	}
-    }
-  delete photonVector;
   if (localEnergyDeposit < 0)
     {
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversion.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversion.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversion.cc	(revision 1055)
@@ -97,4 +97,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeGammaConversion is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
  
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversionModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversionModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeGammaConversionModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeGammaConversionModel.cc,v 1.2 2008/12/04 14:09:36 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeGammaConversionModel.cc,v 1.4 2009/05/19 14:57:01 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -32,4 +32,10 @@
 // --------
 // 06 Oct 2008   L Pandola    Migration from process to model 
+// 17 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - do not apply production threshold on level of the model
+// 19 May 2009   L Pandola    Explicitely set to zero pointers deleted in 
+//                            Initialise(), since they might be checked later on
 //
 
@@ -56,5 +62,5 @@
   fSmallEnergy = 1.1*MeV;
 
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
   //
@@ -89,23 +95,10 @@
       crossSectionHandler->Clear();
       delete crossSectionHandler;
+      crossSectionHandler = 0;
     }
   
-  //Check energy limits
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopeGammaConversionModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopeGammaConversionModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
-
   //Re-initialize cross section handler
   crossSectionHandler = new G4CrossSectionHandler();
-  crossSectionHandler->Initialise(0,LowEnergyLimit(),HighEnergyLimit(),400);
+  crossSectionHandler->Initialise(0,fIntrinsicLowEnergyLimit,HighEnergyLimit(),400);
   crossSectionHandler->Clear();
   G4String crossSectionFile = "penelope/pp-cs-pen-";
@@ -117,16 +110,14 @@
     G4cout << "Loaded cross section files for PenelopeGammaConversion" << G4endl;
 
-  G4cout << "Penelope Gamma Conversion model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / MeV << " MeV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if (verboseLevel > 0) { 
+    G4cout << "Penelope Gamma Conversion model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / MeV << " MeV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -151,10 +142,10 @@
 
   G4int iZ = (G4int) Z;
-  if (!crossSectionHandler)
-    {
-      G4cout << "G4PenelopeGammaConversionModel::ComputeCrossSectionPerAtom" << G4endl;
-      G4cout << "The cross section handler is not correctly initialized" << G4endl;
-      G4Exception();
-    }
+  //  if (!crossSectionHandler) //VI: should not be checked in run time
+  //  {
+  //    G4cout << "G4PenelopeGammaConversionModel::ComputeCrossSectionPerAtom" << G4endl;
+  //    G4cout << "The cross section handler is not correctly initialized" << G4endl;
+  //    G4Exception();
+  //  }
   G4double cs = crossSectionHandler->FindValue(iZ,energy);
 
@@ -167,9 +158,10 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-void G4PenelopeGammaConversionModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
-					      const G4MaterialCutsCouple* couple,
-					      const G4DynamicParticle* aDynamicGamma,
-					      G4double,
-					      G4double)
+void 
+G4PenelopeGammaConversionModel::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
+						  const G4MaterialCutsCouple* couple,
+						  const G4DynamicParticle* aDynamicGamma,
+						  G4double,
+						  G4double)
 {
   //
@@ -191,11 +183,13 @@
   G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
 
-  if (photonEnergy <= LowEnergyLimit())
-  {
-      fParticleChange->ProposeTrackStatus(fStopAndKill);
-      fParticleChange->SetProposedKineticEnergy(0.);
+  // Always kill primary
+  fParticleChange->ProposeTrackStatus(fStopAndKill);
+  fParticleChange->SetProposedKineticEnergy(0.);
+
+  if (photonEnergy <= fIntrinsicLowEnergyLimit)
+    {
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
       return ;
-  }
+    }
 
   G4ParticleMomentum photonDirection = aDynamicGamma->GetMomentumDirection();
@@ -310,11 +304,7 @@
 
   //Generate explicitely the electron in the pair, only if it is > threshold
-  const G4ProductionCutsTable* theCoupleTable=
-    G4ProductionCutsTable::GetProductionCutsTable();
-  size_t indx = couple->GetIndex();
-  G4double cutE = (*(theCoupleTable->GetEnergyCutsVector(1)))[indx];
-  //G4double cutP = (*(theCoupleTable->GetEnergyCutsVector(2)))[indx];
-
-  if (electronKineEnergy > cutE)
+  //VI: applying cut here provides inconsistency 
+
+  if (electronKineEnergy > 0.0)
     {
       G4ThreeVector electronDirection ( dirX_el, dirY_el, dirZ_el);
@@ -333,5 +323,6 @@
   //Generate the positron. Real particle in any case, because it will annihilate. If below
   //threshold, produce it at rest
-  if (positronKineEnergy < cutE)
+  // VI: here there was a bug - positron and electron cuts are different
+  if (positronKineEnergy < 0.0)
     {
       localEnergyDeposit += positronKineEnergy;
@@ -344,11 +335,8 @@
   fvect->push_back(positron);
 
-  //Update the status of the primary gamma (kill it)
-  fParticleChange->SetProposedKineticEnergy(0.);
+  //Add rest of energy to the local energy deposit
   fParticleChange->ProposeLocalEnergyDeposit(localEnergyDeposit);
-  fParticleChange->ProposeTrackStatus(fStopAndKill);
-
   
- if (verboseLevel > 1)
+  if (verboseLevel > 1)
     {
       G4cout << "-----------------------------------------------------------" << G4endl;
@@ -357,5 +345,6 @@
       G4cout << "-----------------------------------------------------------" << G4endl;
       if (electronKineEnergy)
-	G4cout << "Electron (explicitely produced) " << electronKineEnergy/keV << " keV" << G4endl;
+	G4cout << "Electron (explicitely produced) " << electronKineEnergy/keV << " keV" 
+	       << G4endl;
       if (positronKineEnergy)
 	G4cout << "Positron (not at rest) " << positronKineEnergy/keV << " keV" << G4endl;
@@ -373,10 +362,9 @@
 				      localEnergyDeposit+2.0*electron_mass_c2-photonEnergy);
       if (energyDiff > 0.05*keV)
-	G4cout << "Warning from G4PenelopeGammaConversion: problem with energy conservation: " << 
-	  (electronKineEnergy+positronKineEnergy+
-	   localEnergyDeposit+2.0*electron_mass_c2)/keV << " keV (final) vs. " << 
-	  photonEnergy/keV << " keV (initial)" << G4endl;
-    }
-  
+	G4cout << "Warning from G4PenelopeGammaConversion: problem with energy conservation: " 
+	       << (electronKineEnergy+positronKineEnergy+
+		   localEnergyDeposit+2.0*electron_mass_c2)/keV 
+	       << " keV (final) vs. " << photonEnergy/keV << " keV (initial)" << G4endl;
+    } 
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisation.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisation.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisation.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeIonisation.cc,v 1.19 2006/06/29 19:40:49 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeIonisation.cc,v 1.20 2009/05/02 09:59:16 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 // 
 // --------------------------------------------------------------
@@ -93,4 +93,14 @@
   ReadData(); //Read data from file
  
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeIonisation is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
+
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisationModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisationModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeIonisationModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeIonisationModel.cc,v 1.2 2008/12/05 09:15:43 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeIonisationModel.cc,v 1.5 2009/05/19 14:57:01 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -32,4 +32,11 @@
 // --------
 // 26 Nov 2008   L Pandola    Migration from process to model 
+// 17 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - added MinEnergyCut method
+//                  - do not change track status
+// 19 May 2009   L Pandola    Explicitely set to zero pointers deleted in 
+//                            Initialise(), since they might be checked later on
 //
 
@@ -48,5 +55,4 @@
 #include "G4CrossSectionHandler.hh"
 #include "G4AtomicDeexcitation.hh"
-#include "G4ProcessManager.hh"
 #include "G4VEMDataSet.hh"
 
@@ -68,5 +74,5 @@
   fIntrinsicLowEnergyLimit = 100.0*eV;
   fIntrinsicHighEnergyLimit = 100.0*GeV;
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
   //
@@ -136,6 +142,4 @@
     G4cout << "Calling G4PenelopeIonisationModel::Initialise()" << G4endl;
 
-  InitialiseElementSelectors(particle,cuts);
-
   //Delete and re-initialize the cross section handler
   if (crossSectionHandler)
@@ -143,4 +147,5 @@
       crossSectionHandler->Clear();
       delete crossSectionHandler;
+      crossSectionHandler = 0;
     }
 
@@ -148,21 +153,11 @@
     {
       for (size_t i=0; i<theXSTable->size(); i++)
-	delete (*theXSTable)[i];
+	{
+	  delete (*theXSTable)[i];
+	  (*theXSTable)[i] = 0;
+	}
       delete theXSTable;
+      theXSTable = 0;
      }
-
-
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopeIonisationModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopeIonisationModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
 
   crossSectionHandler = new G4CrossSectionHandler();
@@ -177,21 +172,20 @@
   //This is used to retrieve cross section values later on
   crossSectionHandler->BuildMeanFreePathForMaterials();
-  
- if (verboseLevel > 2) 
+
+  InitialiseElementSelectors(particle,cuts);
+  
+  if (verboseLevel > 2) 
     G4cout << "Loaded cross section files for PenelopeIonisationModel" << G4endl;
 
-
-  G4cout << "Penelope Ionisation model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if (verboseLevel > 2) {
+    G4cout << "Penelope Ionisation model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / keV << " keV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForLoss*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForLoss();
+  fParticleChange = GetParticleChangeForLoss();
   isInitialised = true; 
 }
@@ -227,4 +221,6 @@
   SetupForMaterial(theParticle, material, energy);
 
+  // VI - should be check at initialisation not in run time
+  /*
   if (!crossSectionHandler)
     {
@@ -233,5 +229,5 @@
       G4Exception();
     }
-  
+  */  
   if (!theXSTable)
     {
@@ -244,5 +240,5 @@
       theXSTable = BuildCrossSectionTable(theParticle); 
     }
-
+  
   G4double totalCross = 0.0; 
   G4double cross = 0.0;
@@ -356,6 +352,5 @@
 					      const G4MaterialCutsCouple* couple,
 					      const G4DynamicParticle* aDynamicParticle,
-					      G4double,
-					      G4double)
+					      G4double cutE, G4double)
 {
   // Penelope model to sample the final state following an hard inelastic interaction.
@@ -397,14 +392,8 @@
   const G4ParticleDefinition* theParticle = aDynamicParticle->GetDefinition();
 
-  if (kineticEnergy0 <= LowEnergyLimit())
+  if (kineticEnergy0 <= fIntrinsicLowEnergyLimit)
   {
     fParticleChange->SetProposedKineticEnergy(0.);
     fParticleChange->ProposeLocalEnergyDeposit(kineticEnergy0);
-    //Check if there are AtRest processes
-    if (theParticle->GetProcessManager()->GetAtRestProcessVector()->size()) 
-      //In this case there is at least one AtRest process
-      fParticleChange->ProposeTrackStatus(fStopButAlive);
-    else
-      fParticleChange->ProposeTrackStatus(fStopAndKill);
     return ;
   }
@@ -430,15 +419,8 @@
   G4int iZ = SampleRandomAtom(couple,kineticEnergy0);
   
-  G4double cutForLowEnergySecondaryParticles = 250.0*eV;
   const G4ProductionCutsTable* theCoupleTable=
     G4ProductionCutsTable::GetProductionCutsTable();
   size_t indx = couple->GetIndex();
-  G4double cutE = (*(theCoupleTable->GetEnergyCutsVector(1)))[indx];
   G4double cutG = (*(theCoupleTable->GetEnergyCutsVector(0)))[indx];
-
-  //Production cut for delta-rays (electrons)
-  cutE = std::max(cutForLowEnergySecondaryParticles,cutE);
-  //Production cut for gamma (fluorescence)
-  cutG = std::max(cutForLowEnergySecondaryParticles,cutG);
   
   if (verboseLevel > 2)
@@ -474,11 +456,5 @@
   else
     {
-      fParticleChange->ProposeMomentumDirection(electronDirection1);
-      fParticleChange->SetProposedKineticEnergy(0.*eV);
-      if (theParticle->GetProcessManager()->GetAtRestProcessVector()->size()) 
-	//In this case there is at least one AtRest process
-	fParticleChange->ProposeTrackStatus(fStopButAlive);
-      else
-	fParticleChange->ProposeTrackStatus(fStopAndKill);
+      fParticleChange->SetProposedKineticEnergy(0.);
     }
   
@@ -517,36 +493,40 @@
     //Penelope subtracted the fluorescence, but one has to match the databases
     eKineticEnergy = energySecondary+ioniEnergy-bindingEnergy;
-    
-  //VERIFICA QUI LA STORIA DEL LOCAL ENERGY DEPOSIT!
-
+  
   G4double localEnergyDeposit = ionEnergy; 
   G4double energyInFluorescence = 0.0*eV;
 
-   std::vector<G4DynamicParticle*> *photonVector = 0;
-  if (fUseAtomicDeexcitation)
-    {
-      if (iZ>5 && (ionEnergy > cutG || ionEnergy > cutE))
+  if(DeexcitationFlag() && iZ > 5) 
+    {
+      if (ionEnergy > cutG || ionEnergy > cutE)
 	{
-	  photonVector = deexcitationManager.GenerateParticles(iZ,shellId);
-	  //Check for single photons if they are above threshold
-	  for (size_t k=0;k<photonVector->size();k++)
+	  deexcitationManager.SetCutForSecondaryPhotons(cutG);
+	  deexcitationManager.SetCutForAugerElectrons(cutE);
+	  std::vector<G4DynamicParticle*> *photonVector =
+	    deexcitationManager.GenerateParticles(iZ,shellId);
+	  //Check for secondaries
+          if(photonVector) 
 	    {
-	      G4DynamicParticle* aPhoton = (*photonVector)[k];
-	      if (aPhoton)
+	      for (size_t k=0;k<photonVector->size();k++)
 		{
-		  G4double itsCut = cutG;
-		  if (aPhoton->GetDefinition() == G4Electron::Electron()) itsCut = cutE;
-		  G4double itsEnergy = aPhoton->GetKineticEnergy();
-		  if (itsEnergy > itsCut && itsEnergy <= ionEnergy)
+		  G4DynamicParticle* aPhoton = (*photonVector)[k];
+		  if (aPhoton)
 		    {
-		      localEnergyDeposit -= itsEnergy;
-		      energyInFluorescence += itsEnergy;
-		    }
-		  else
-		    {
-		      delete aPhoton;
-		      (*photonVector)[k] = 0;
+		      G4double itsEnergy = aPhoton->GetKineticEnergy();
+		      if (itsEnergy <= localEnergyDeposit)
+			{
+			  if(aPhoton->GetDefinition() == G4Gamma::Gamma())
+			    energyInFluorescence += itsEnergy;
+			  localEnergyDeposit -= itsEnergy;
+			  fvect->push_back(aPhoton);
+			}
+		      else
+			{
+			  delete aPhoton;
+			  (*photonVector)[k] = 0;
+			}
 		    }
 		}
+	      delete photonVector;
 	    }
 	}
@@ -567,18 +547,4 @@
   fvect->push_back(deltaElectron);
   
-  //Generate fluorescence, if it is the case
-  //This block is executed only if there is at least one secondary photon produced by 
-  //G4AtomicDeexcitation
-  if (photonVector)
-    {
-      for (size_t ll=0;ll<photonVector->size();ll++)
-	if ((*photonVector)[ll])
-	  {
-	    G4DynamicParticle* aFluorescencePhoton = (*photonVector)[ll];
-	    fvect->push_back(aFluorescencePhoton);
-	  }
-    }
-  delete photonVector;
-
   if (localEnergyDeposit < 0)
     {
@@ -753,6 +719,9 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-void G4PenelopeIonisationModel::CalculateDiscreteForElectrons(G4double kinEnergy,G4double cutoffEnergy,
-							      G4int Z,G4double electronVolumeDensity)
+void 
+G4PenelopeIonisationModel::CalculateDiscreteForElectrons(G4double kinEnergy,
+							 G4double cutoffEnergy,
+							 G4int Z,
+							 G4double electronVolumeDensity)
 {
   if (verboseLevel > 2)
@@ -1030,6 +999,9 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
- void G4PenelopeIonisationModel::CalculateDiscreteForPositrons(G4double kinEnergy,G4double cutoffEnergy,
-							       G4int Z,G4double electronVolumeDensity)
+ void 
+G4PenelopeIonisationModel::CalculateDiscreteForPositrons(G4double kinEnergy,
+							 G4double cutoffEnergy,
+							 G4int Z,
+							 G4double electronVolumeDensity)
 {
   kineticEnergy1=kinEnergy;
@@ -1287,8 +1259,9 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-G4double G4PenelopeIonisationModel::CalculateCrossSectionsRatio(G4double kinEnergy,
-								G4double cutoffEnergy,
-								G4int Z, G4double electronVolumeDensity,
-								const G4ParticleDefinition* theParticle)
+G4double 
+G4PenelopeIonisationModel::CalculateCrossSectionsRatio(G4double kinEnergy,
+						       G4double cutoffEnergy,
+						       G4int Z, G4double electronVolumeDensity,
+						       const G4ParticleDefinition* theParticle)
 {
   //Constants
@@ -1320,8 +1293,9 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-std::pair<G4double,G4double> G4PenelopeIonisationModel::CrossSectionsRatioForElectrons(G4double kineticEnergy,
-										       G4double resEnergy,
-										       G4double densityCorrection,
-										       G4double cutoffEnergy)
+std::pair<G4double,G4double> 
+G4PenelopeIonisationModel::CrossSectionsRatioForElectrons(G4double kineticEnergy,
+							  G4double resEnergy,
+							  G4double densityCorrection,
+							  G4double cutoffEnergy)
 {
   std::pair<G4double,G4double> theResult(0.,0.);
@@ -1410,8 +1384,9 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-std::pair<G4double,G4double> G4PenelopeIonisationModel::CrossSectionsRatioForPositrons(G4double kineticEnergy,
-										       G4double resEnergy,
-										       G4double densityCorrection,
-										       G4double cutoffEnergy)
+std::pair<G4double,G4double> 
+G4PenelopeIonisationModel::CrossSectionsRatioForPositrons(G4double kineticEnergy,
+							  G4double resEnergy,
+							  G4double densityCorrection,
+							  G4double cutoffEnergy)
 {
 
@@ -1640,6 +1615,6 @@
 #include "G4EMDataSet.hh"
 
-std::vector<G4VEMDataSet*>* G4PenelopeIonisationModel::BuildCrossSectionTable(const 
-									      G4ParticleDefinition* theParticle)
+std::vector<G4VEMDataSet*>* 
+G4PenelopeIonisationModel::BuildCrossSectionTable(const G4ParticleDefinition* theParticle)
 {
   std::vector<G4VEMDataSet*>* set = new std::vector<G4VEMDataSet*>;
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectric.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectric.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectric.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4PenelopePhotoElectric.cc,v 1.12 2006/06/29 19:40:51 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopePhotoElectric.cc,v 1.14 2009/05/02 09:59:17 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: L. Pandola
@@ -37,4 +37,9 @@
 // 31 May 2005  L. Pandola  Added Sauter formula for the sampling of 
 //                          the electron direction
+// 08 Jan 2009  L. Pandola  Check shell index to avoid mismatch between 
+//                          the Penelope cross section database and the 
+//                          G4AtomicTransitionManager database. It suppresses 
+//                          a warning from G4AtomicTransitionManager only. 
+//                          Results are unchanged.
 // --------------------------------------------------------------
 
@@ -89,4 +94,15 @@
 	     << G4endl;
     }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopePhotoElectric is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
+
 }
 
@@ -146,4 +162,15 @@
   // Retrieve the corresponding identifier and binding energy of the selected shell
   const G4AtomicTransitionManager* transitionManager = G4AtomicTransitionManager::Instance();
+
+  //The number of shell cross section possibly reported in the Penelope database 
+  //might be different from the number of shells in the G4AtomicTransitionManager
+  //(namely, Penelope may contain more shell, especially for very light elements).
+  //In order to avoid a warning message from the G4AtomicTransitionManager, I 
+  //add this protection. Results are anyway changed, because when G4AtomicTransitionManager
+  //has a shellID>maxID, it sets the shellID to the last valid shell. 
+  size_t numberOfShells = (size_t) transitionManager->NumberOfShells(Z);
+  if (shellIndex >= numberOfShells)
+    shellIndex = numberOfShells-1;
+
   const G4AtomicShell* shell = transitionManager->Shell(Z,shellIndex);
   G4double bindingEnergy = shell->BindingEnergy();
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectricModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectricModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopePhotoElectricModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopePhotoElectricModel.cc,v 1.2 2008/12/04 14:09:36 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopePhotoElectricModel.cc,v 1.6 2009/05/19 14:57:01 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -31,5 +31,17 @@
 // History:
 // --------
-// 08 Oct 2008   L Pandola    Migration from process to model 
+// 08 Oct 2008   L Pandola  Migration from process to model 
+// 08 Jan 2009   L Pandola  Check shell index to avoid mismatch between 
+//                          the Penelope cross section database and the 
+//                          G4AtomicTransitionManager database. It suppresses 
+//                          a warning from G4AtomicTransitionManager only. 
+//                          Results are unchanged.
+// 25 Mar 2009   L Pandola  Small fix to avoid wrong energy-violation warnings
+// 17 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+//                  - do not apply production threshold on secondaries
+// 19 May 2009   L Pandola    Explicitely set to zero pointers deleted in 
+//                            Initialise(), since they might be checked later on
 //
 
@@ -58,5 +70,5 @@
   fIntrinsicLowEnergyLimit = 100.0*eV;
   fIntrinsicHighEnergyLimit = 100.0*GeV;
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
   //
@@ -90,4 +102,5 @@
       crossSectionHandler->Clear();
       delete crossSectionHandler;
+      crossSectionHandler = 0;
     }
   if (shellCrossSectionHandler)
@@ -95,22 +108,6 @@
       shellCrossSectionHandler->Clear();
       delete shellCrossSectionHandler;
-    }
-
-  //Check energy limits
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopePhotoElectricModel: low energy limit increased from " <<
-        LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << 
-	G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopePhotoElectricModel: high energy limit decreased from " <<
-        HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" 
-	     << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
-
+      shellCrossSectionHandler =0;
+    }
 
   //Re-initialize cross section handlers
@@ -129,16 +126,14 @@
     G4cout << "Loaded cross section files for PenelopePhotoElectric" << G4endl;
 
-  G4cout << "Penelope Photo-Electric model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / MeV << " MeV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+  if (verboseLevel > 0) { 
+    G4cout << "Penelope Photo-Electric model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / MeV << " MeV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -162,10 +157,10 @@
 
   G4int iZ = (G4int) Z;
-  if (!crossSectionHandler)
-    {
-      G4cout << "G4PenelopePhotoElectricModel::ComputeCrossSectionPerAtom" << G4endl;
-      G4cout << "The cross section handler is not correctly initialized" << G4endl;
-      G4Exception();
-    }
+  //  if (!crossSectionHandler) // VI: should not be 
+  //  {
+  //    G4cout << "G4PenelopePhotoElectricModel::ComputeCrossSectionPerAtom" << G4endl;
+  //    G4cout << "The cross section handler is not correctly initialized" << G4endl;
+  //    G4Exception();
+  //  }
   G4double cs = crossSectionHandler->FindValue(iZ,energy);
  
@@ -203,11 +198,13 @@
   G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
 
-  if (photonEnergy <= LowEnergyLimit())
-  {
-      fParticleChange->ProposeTrackStatus(fStopAndKill);
-      fParticleChange->SetProposedKineticEnergy(0.);
+  // always kill primary
+  fParticleChange->ProposeTrackStatus(fStopAndKill);
+  fParticleChange->SetProposedKineticEnergy(0.);
+
+  if (photonEnergy <= fIntrinsicLowEnergyLimit)
+    {
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
       return ;
-  }
+    }
 
   G4ParticleMomentum photonDirection = aDynamicGamma->GetMomentumDirection();
@@ -217,5 +214,6 @@
     G4cout << "Going to select element in " << couple->GetMaterial()->GetName() << G4endl;
   //use crossSectionHandler instead of G4EmElementSelector because in this case 
-  //the dimension of the table is equal to the dimension of the database (less interpolation errors)
+  //the dimension of the table is equal to the dimension of the database 
+  //(less interpolation errors)
   G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
   if (verboseLevel > 2)
@@ -227,4 +225,15 @@
   // Retrieve the corresponding identifier and binding energy of the selected shell
   const G4AtomicTransitionManager* transitionManager = G4AtomicTransitionManager::Instance();
+
+  //The number of shell cross section possibly reported in the Penelope database 
+  //might be different from the number of shells in the G4AtomicTransitionManager
+  //(namely, Penelope may contain more shell, especially for very light elements).
+  //In order to avoid a warning message from the G4AtomicTransitionManager, I 
+  //add this protection. Results are anyway changed, because when G4AtomicTransitionManager
+  //has a shellID>maxID, it sets the shellID to the last valid shell. 
+  size_t numberOfShells = (size_t) transitionManager->NumberOfShells(Z);
+  if (shellIndex >= numberOfShells)
+    shellIndex = numberOfShells-1;
+
   const G4AtomicShell* shell = transitionManager->Shell(Z,shellIndex);
   G4double bindingEnergy = shell->BindingEnergy();
@@ -236,11 +245,4 @@
   G4double eKineticEnergy = photonEnergy - bindingEnergy;
   
-  G4double cutForLowEnergySecondaryParticles = 250.0*eV;
-  const G4ProductionCutsTable* theCoupleTable=
-    G4ProductionCutsTable::GetProductionCutsTable();
-  size_t indx = couple->GetIndex();
-  G4double cutE = (*(theCoupleTable->GetEnergyCutsVector(1)))[indx];
-  cutE = std::max(cutForLowEnergySecondaryParticles,cutE);
-
   // There may be cases where the binding energy of the selected shell is > photon energy
   // In such cases do not generate secondaries
@@ -248,81 +250,75 @@
     {
       //Now check if the electron is above cuts: if so, it is created explicitely
-      if (eKineticEnergy > cutE)
-	{
-	  // The electron is created
-	  // Direction sampled from the Sauter distribution
-	  G4double cosTheta = SampleElectronDirection(eKineticEnergy);
-	  G4double sinTheta = std::sqrt(1-cosTheta*cosTheta);
-	  G4double phi = twopi * G4UniformRand() ;
-	  G4double dirx = sinTheta * std::cos(phi);
-	  G4double diry = sinTheta * std::sin(phi);
-	  G4double dirz = cosTheta ;
-	  G4ThreeVector electronDirection(dirx,diry,dirz); //electron direction
-	  electronDirection.rotateUz(photonDirection);
-	  G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(), 
-							       electronDirection, 
-							       eKineticEnergy);
-	  fvect->push_back(electron);
-	} 
-      else 
-	localEnergyDeposit += eKineticEnergy;    
-    }
+      //VI: checking cut here provides inconsistency in testing
+      //      if (eKineticEnergy > cutE)
+      // The electron is created
+      // Direction sampled from the Sauter distribution
+      G4double cosTheta = SampleElectronDirection(eKineticEnergy);
+      G4double sinTheta = std::sqrt(1-cosTheta*cosTheta);
+      G4double phi = twopi * G4UniformRand() ;
+      G4double dirx = sinTheta * std::cos(phi);
+      G4double diry = sinTheta * std::sin(phi);
+      G4double dirz = cosTheta ;
+      G4ThreeVector electronDirection(dirx,diry,dirz); //electron direction
+      electronDirection.rotateUz(photonDirection);
+      G4DynamicParticle* electron = new G4DynamicParticle (G4Electron::Electron(), 
+							   electronDirection, 
+							   eKineticEnergy);
+      fvect->push_back(electron);
+    } 
+  //  else 
+  //  {
+  //    localEnergyDeposit += eKineticEnergy;    
+  //    eKineticEnergy = 0;
+  //  }
+  //  }
   else
+    {
       bindingEnergy = photonEnergy;
-    
-  
+    }
   G4double energyInFluorescence = 0; //testing purposes
 
   //Now, take care of fluorescence, if required
-  if (fUseAtomicDeexcitation)
-    {
+  if(DeexcitationFlag() && Z > 5) 
+    {
+      const G4ProductionCutsTable* theCoupleTable=
+	G4ProductionCutsTable::GetProductionCutsTable();
+      size_t indx = couple->GetIndex();
       G4double cutG = (*(theCoupleTable->GetEnergyCutsVector(0)))[indx];
-      cutG = std::min(cutForLowEnergySecondaryParticles,cutG);
+      G4double cutE = (*(theCoupleTable->GetEnergyCutsVector(1)))[indx];
       
-      std::vector<G4DynamicParticle*>* photonVector = 0;
-
       // Protection to avoid generating photons in the unphysical case of 
       // shell binding energy > photon energy
-      if (Z > 5  && (bindingEnergy > cutG || bindingEnergy > cutE))
+      if (bindingEnergy > cutG || bindingEnergy > cutE)
 	{
-	  photonVector = deexcitationManager.GenerateParticles(Z,shellId); 
-	  //Check for single photons (if they are above threshold)
-	  for (size_t k=0; k< photonVector->size(); k++)
+	  deexcitationManager.SetCutForSecondaryPhotons(cutG);
+	  deexcitationManager.SetCutForAugerElectrons(cutE);
+	  std::vector<G4DynamicParticle*>* photonVector = 
+	    deexcitationManager.GenerateParticles(Z,shellId); 
+	  //Check for secondaries
+          if(photonVector) 
 	    {
-	      G4DynamicParticle* aPhoton = (*photonVector)[k];
-	      if (aPhoton)
+	      for (size_t k=0; k< photonVector->size(); k++)
 		{
-		  G4double itsCut = cutG;
-		  if(aPhoton->GetDefinition() == G4Electron::Electron()) itsCut = cutE;
-		  G4double itsEnergy = aPhoton->GetKineticEnergy();
-		  
-		  if (itsEnergy > itsCut && itsEnergy <= bindingEnergy)
+		  G4DynamicParticle* aPhoton = (*photonVector)[k];
+		  if (aPhoton)
 		    {
-		      // Local energy deposit is given as the sum of the 
-		      // energies of incident photons minus the energies
-		      // of the outcoming fluorescence photons
-		      bindingEnergy -= itsEnergy;
-		      
-		    }
-		  else
-		    { 
-		      (*photonVector)[k] = 0;
+		      G4double itsEnergy = aPhoton->GetKineticEnergy();
+		      if (itsEnergy <= bindingEnergy)
+			{
+			  if(aPhoton->GetDefinition() == G4Gamma::Gamma())
+			    energyInFluorescence += itsEnergy;
+			  bindingEnergy -= itsEnergy;
+			  fvect->push_back(aPhoton);
+			}
+		      else
+			{
+			  delete aPhoton;
+			  (*photonVector)[k] = 0;
+			}
 		    }
 		}
+	      delete photonVector;
 	    }
-	}
-      //Register valid secondaries
-      if (photonVector)
-	{
-	  for ( size_t ll = 0; ll <photonVector->size(); ll++) 
-	    {
-	      G4DynamicParticle* aPhoton = (*photonVector)[ll];
-	      if (aPhoton) 
-		{
-		  energyInFluorescence += aPhoton->GetKineticEnergy();
-		  fvect->push_back(aPhoton);
-		}
-	    }
-	  delete photonVector;
 	}
     }
@@ -330,5 +326,4 @@
   localEnergyDeposit += bindingEnergy;
       
-
   if (localEnergyDeposit < 0)
     {
@@ -339,8 +334,5 @@
     }
 
- //Update the status of the primary gamma (kill it)
-  fParticleChange->SetProposedKineticEnergy(0.);
   fParticleChange->ProposeLocalEnergyDeposit(localEnergyDeposit);
-  fParticleChange->ProposeTrackStatus(fStopAndKill);
 
   if (verboseLevel > 1)
@@ -360,8 +352,10 @@
   if (verboseLevel > 0)
     {
-      G4double energyDiff = std::fabs(eKineticEnergy+energyInFluorescence+localEnergyDeposit-photonEnergy);
+      G4double energyDiff = 
+	std::fabs(eKineticEnergy+energyInFluorescence+localEnergyDeposit-photonEnergy);
       if (energyDiff > 0.05*keV)
 	G4cout << "Warning from G4PenelopePhotoElectric: problem with energy conservation: " << 
-	  (eKineticEnergy+energyInFluorescence+localEnergyDeposit)/keV << " keV (final) vs. " << 
+	  (eKineticEnergy+energyInFluorescence+localEnergyDeposit)/keV 
+	       << " keV (final) vs. " << 
 	  photonEnergy/keV << " keV (initial)" << G4endl;
     }
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleigh.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleigh.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleigh.cc	(revision 1055)
@@ -26,6 +26,6 @@
 // --------------------------------------------------------------------
 //
-// $Id: G4PenelopeRayleigh.cc,v 1.15 2007/09/03 09:43:14 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeRayleigh.cc,v 1.17 2009/05/02 09:59:17 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: L. Pandola (luciano.pandola@cern.ch)
@@ -41,4 +41,7 @@
 // 03 Sep 2007 L.Pandola      Bug fix for the filling of physics table for 
 //                            compounds defined by the mass fraction (bug #965)
+// 02 Apr 2009 L.Pandola      Bux fixed in the calculation of mfp for compound 
+//                            materials defined as fraction of mass 
+//                            (reported by Zhang Qiwei)
 // --------------------------------------------------------------------
 
@@ -93,4 +96,14 @@
 	      << G4endl;
      }
+
+   G4cout << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "   The class G4PenelopeRayleigh is NOT SUPPORTED ANYMORE. " << G4endl;
+   G4cout << "   It will be REMOVED with the next major release of Geant4. " << G4endl;
+   G4cout << "   Please consult: https://twiki.cern.ch/twiki/bin/view/Geant4/LoweProcesses" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << "*******************************************************************************" << G4endl;
+   G4cout << G4endl;
 }
 
@@ -162,5 +175,7 @@
 	  const G4double* vector_of_atoms = material->GetVecNbOfAtomsPerVolume();
 	  const G4int* stechiometric = material->GetAtomsVector();
-	  G4double density = vector_of_atoms[iright]; //non-bound molecules (default)
+	  //cs is the cross section _per atom_ in the case of compounds, while it is 
+	  //_per molecule_ in the case of molecules
+	  G4double density = material->GetTotNbOfAtomsPerVolume();  //non-bound molecules (default)
 	  if (stechiometric)
 	    {
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleighModel.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleighModel.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4PenelopeRayleighModel.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PenelopeRayleighModel.cc,v 1.2 2008/12/04 14:09:36 pandola Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4PenelopeRayleighModel.cc,v 1.4 2009/05/19 14:57:01 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Luciano Pandola
@@ -32,4 +32,9 @@
 // --------
 // 14 Oct 2008   L Pandola    Migration from process to model 
+// 17 Apr 2009   V Ivanchenko Cleanup initialisation and generation of secondaries:
+//                  - apply internal high-energy limit only in constructor 
+//                  - do not apply low-energy limit (default is 0)
+// 19 May 2009   L Pandola    Explicitely set to zero pointers deleted in 
+//                            PrepareConstants(), since they might be checked later on
 //
 
@@ -56,5 +61,5 @@
   fIntrinsicLowEnergyLimit = 100.0*eV;
   fIntrinsicHighEnergyLimit = 100.0*GeV;
-  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
+  //  SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
   //
@@ -90,29 +95,15 @@
     G4cout << "Calling G4PenelopeRayleighModel::Initialise()" << G4endl;
 
-  if (LowEnergyLimit() < fIntrinsicLowEnergyLimit)
-    {
-      G4cout << "G4PenelopeRayleighModel: low energy limit increased from " << 
-	LowEnergyLimit()/eV << " eV to " << fIntrinsicLowEnergyLimit/eV << " eV" << G4endl;
-      SetLowEnergyLimit(fIntrinsicLowEnergyLimit);
-    }
-  if (HighEnergyLimit() > fIntrinsicHighEnergyLimit)
-    {
-      G4cout << "G4PenelopeRayleighModel: high energy limit decreased from " << 
-	HighEnergyLimit()/GeV << " GeV to " << fIntrinsicHighEnergyLimit/GeV << " GeV" << G4endl;
-      SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
-    }
-
-  G4cout << "Penelope Rayleigh model is initialized " << G4endl
-         << "Energy range: "
-         << LowEnergyLimit() / keV << " keV - "
-         << HighEnergyLimit() / GeV << " GeV"
-         << G4endl;
+
+  if (verboseLevel > 0) {
+    G4cout << "Penelope Rayleigh model is initialized " << G4endl
+	   << "Energy range: "
+	   << LowEnergyLimit() / keV << " keV - "
+	   << HighEnergyLimit() / GeV << " GeV"
+	   << G4endl;
+  }
 
   if(isInitialised) return;
-
-  if(pParticleChange)
-    fParticleChange = reinterpret_cast<G4ParticleChangeForGamma*>(pParticleChange);
-  else
-    fParticleChange = new G4ParticleChangeForGamma();
+  fParticleChange = GetParticleChangeForGamma();
   isInitialised = true;
 }
@@ -120,9 +111,10 @@
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 
-G4double G4PenelopeRayleighModel::CrossSectionPerVolume(const G4Material* material,
-                                           const G4ParticleDefinition* p,
-                                           G4double ekin,
-                                           G4double,
-                                           G4double)
+G4double 
+G4PenelopeRayleighModel::CrossSectionPerVolume(const G4Material* material,
+					       const G4ParticleDefinition* p,
+					       G4double ekin,
+					       G4double,
+					       G4double)
 {
   // Penelope model to calculate the Rayleigh scattering inverse mean 
@@ -233,8 +225,8 @@
 
 void G4PenelopeRayleighModel::SampleSecondaries(std::vector<G4DynamicParticle*>* ,
-					      const G4MaterialCutsCouple* couple,
-					      const G4DynamicParticle* aDynamicGamma,
-					      G4double,
-					      G4double)
+						const G4MaterialCutsCouple* couple,
+						const G4DynamicParticle* aDynamicGamma,
+						G4double,
+						G4double)
 {
   // Penelope model to sample the Rayleigh scattering final state.
@@ -249,14 +241,14 @@
   if (verboseLevel > 3)
     G4cout << "Calling SamplingSecondaries() of G4PenelopeRayleighModel" << G4endl;
-                                                                                                    
+
   G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
  
-  if (photonEnergy0 <= LowEnergyLimit())
-  {
+  if (photonEnergy0 <= fIntrinsicLowEnergyLimit)
+    {
       fParticleChange->ProposeTrackStatus(fStopAndKill);
       fParticleChange->SetProposedKineticEnergy(0.);
       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
       return ;
-  }
+    }
 
   G4ParticleMomentum photonDirection0 = aDynamicGamma->GetMomentumDirection();
@@ -534,8 +526,14 @@
   Xhigh=1e06;
   if (samplingFunction_x)
-    delete samplingFunction_x;
+    {
+      delete samplingFunction_x;
+      samplingFunction_x = 0;
+    }
   if (samplingFunction_xNoLog)
-    delete samplingFunction_xNoLog;
- 
+    {
+      delete samplingFunction_xNoLog;
+      samplingFunction_xNoLog = 0;
+    }
+
   samplingFunction_x = new G4DataVector();
   samplingFunction_xNoLog = new G4DataVector();
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4ShellVacancy.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4ShellVacancy.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4ShellVacancy.cc	(revision 1055)
@@ -50,4 +50,5 @@
 
 {
+/*
   G4int size = xsis.size();
   for (G4int k =0; k<size; k++)
@@ -56,10 +57,10 @@
       xsis[k] = 0;
     }
+*/
 }
 
-void G4ShellVacancy::AddXsiTable(G4VEMDataSet* set)
-
+void G4ShellVacancy::AddXsiTable(G4VEMDataSet* p)
 {
-  xsis.push_back(set);
+  xsis.push_back(p);
 }
 
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationCrossSectionHandler.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationCrossSectionHandler.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationCrossSectionHandler.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4eIonisationCrossSectionHandler.cc,v 1.11 2006/06/29 19:42:00 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4eIonisationCrossSectionHandler.cc,v 1.12 2009/01/29 08:13:34 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // -------------------------------------------------------------------
@@ -42,4 +42,6 @@
 // 19 Jul 2002   VI          Create composite data set for material
 // 21 Jan 2003  V.Ivanchenko Cut per region
+// 28 Jan 2009  L.Pandola    Added public method to make a easier migration of 
+//                           G4LowEnergyIonisation to G4LivermoreIonisationModel
 //
 // -------------------------------------------------------------------
@@ -98,9 +100,10 @@
     G4int nElements = material->GetNumberOfElements();
 
-    if(verbose > 0) {
-      G4cout << "eIonisation CS for " << m << "th material "
-             << material->GetName()
-             << "  eEl= " << nElements << G4endl;
-    }
+    if(verbose > 0) 
+      {
+	G4cout << "eIonisation CS for " << m << "th material "
+	       << material->GetName()
+	       << "  eEl= " << nElements << G4endl;
+      }
 
     G4double tcut  = (*energyCuts)[m];
@@ -129,5 +132,6 @@
             value += cross * p * density;
 
-            if(verbose>0 && m == 0 && e>=1. && e<=0.) {
+	    if(verbose>0 && m == 0 && e>=1. && e<=0.) 
+	    {
               G4cout << "G4eIonCrossSH: e(MeV)= " << e/MeV
                      << " n= " << n
@@ -139,5 +143,5 @@
                      << " Z= " << Z
                      << G4endl;
-	    }
+	      }
 
 	  }
@@ -155,3 +159,18 @@
 }
 
-
+G4double G4eIonisationCrossSectionHandler::GetCrossSectionAboveThresholdForElement(G4double energy,
+										   G4double cutEnergy,
+										   G4int Z)
+{
+  G4int nShells = NumberOfComponents(Z);
+  G4double value = 0.;
+  if(energy > cutEnergy) 
+    {
+      for (G4int n=0; n<nShells; n++) {
+	G4double cross = FindValue(Z, energy, n);
+	G4double p = theParam->Probability(Z, cutEnergy, energy, energy, n);
+	value += cross * p;
+      }
+    }
+  return value;
+}
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationParameters.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationParameters.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationParameters.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4eIonisationParameters.cc,v 1.23 2006/06/29 19:42:02 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4eIonisationParameters.cc,v 1.24 2009/03/23 09:13:44 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // Author: Maria Grazia Pia (Maria.Grazia.Pia@cern.ch)
@@ -41,4 +41,5 @@
 //                          chenged to "ion-..."
 // 17.02.04 V.Ivanchenko    Increase buffer size
+// 23.03.09 L.Pandola       Updated warning message
 //
 // -------------------------------------------------------------------
@@ -199,5 +200,5 @@
     if (! (lsdp->is_open()) ) {
       G4String excep = G4String("G4IonisationParameters - data file: ")
-      + name + G4String(" not found. The version 1.# of G4LEDATA should be used");
+      + name + G4String(" not found. Please check and/or update G4LEDATA");
       G4Exception(excep);
     }
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationSpectrum.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationSpectrum.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4eIonisationSpectrum.cc	(revision 1055)
@@ -24,6 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4eIonisationSpectrum.cc,v 1.25 2006/06/29 19:42:04 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4eIonisationSpectrum.cc,v 1.26 2009/03/23 09:13:44 pandola Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // -------------------------------------------------------------------
@@ -46,4 +46,6 @@
 // 30.05.2002 VI            Update to 24-parameters data
 // 11.07.2002 VI            Fix in integration over spectrum
+// 23.03.2009 LP            Added protection against division by zero (for 
+//                          faulty database files), Bug Report 1042
 //
 // -------------------------------------------------------------------
@@ -125,5 +127,14 @@
   p.push_back((2.0*g - 1.0)/(g*g));
   
-  p[iMax-1] = Function(p[3], p);
+  //Add protection against division by zero: actually in Function() 
+  //parameter p[3] appears in the denominator. Bug report 1042
+  if (p[3] > 0)
+    p[iMax-1] = Function(p[3], p);
+  else
+    {
+      G4cout << "WARNING: G4eIonisationSpectrum::Probability "
+	     << "parameter p[3] <= 0. G4LEDATA dabatase might be corrupted for Z = " 
+	     << Z << ". Please check and/or update it " << G4endl;      
+    }
 
   if(e >= 1. && e <= 0. && Z == 4) p.push_back(0.0);
@@ -210,5 +221,15 @@
   p.push_back((2.0*g - 1.0)/(g*g));
 
-  p[iMax-1] = Function(p[3], p);
+  
+  //Add protection against division by zero: actually in Function() 
+  //parameter p[3] appears in the denominator. Bug report 1042
+  if (p[3] > 0)
+    p[iMax-1] = Function(p[3], p);
+  else
+    {
+      G4cout << "WARNING: G4eIonisationSpectrum::AverageEnergy "
+	     << "parameter p[3] <= 0. G4LEDATA dabatase might be corrupted for Z = " 
+	     << Z << ". Please check and/or update it " << G4endl;      
+    }
     
   G4double val = AverageValue(x1, x2, p);
@@ -288,5 +309,15 @@
   p.push_back((2.0*g - 1.0)/(g*g));
 
-  p[iMax-1] = Function(p[3], p);
+  
+  //Add protection against division by zero: actually in Function() 
+  //parameter p[3] appears in the denominator. Bug report 1042
+  if (p[3] > 0)
+    p[iMax-1] = Function(p[3], p);
+  else
+    {
+      G4cout << "WARNING: G4eIonisationSpectrum::SampleSpectrum "
+	     << "parameter p[3] <= 0. G4LEDATA dabatase might be corrupted for Z = " 
+	     << Z << ". Please check and/or update it " << G4endl;      
+    }
 
   G4double aria1 = 0.0;
Index: trunk/source/processes/electromagnetic/lowenergy/src/G4hLowEnergyLoss.cc
===================================================================
--- trunk/source/processes/electromagnetic/lowenergy/src/G4hLowEnergyLoss.cc	(revision 1013)
+++ trunk/source/processes/electromagnetic/lowenergy/src/G4hLowEnergyLoss.cc	(revision 1055)
@@ -25,6 +25,6 @@
 //
 //
-// $Id: G4hLowEnergyLoss.cc,v 1.27 2008/06/20 19:54:03 sincerti Exp $
-// GEANT4 tag $Name: geant4-09-02 $
+// $Id: G4hLowEnergyLoss.cc,v 1.28 2009/02/20 10:49:54 sincerti Exp $
+// GEANT4 tag $Name: geant4-09-03-beta-cand-01 $
 //
 // -----------------------------------------------------------
@@ -145,5 +145,6 @@
 G4double G4hLowEnergyLoss::HighestKineticEnergy= 100.*GeV;
 G4int    G4hLowEnergyLoss::TotBin = 360;
-G4double G4hLowEnergyLoss::RTable,G4hLowEnergyLoss::LOGRTable;
+G4double G4hLowEnergyLoss::RTable =1.1;
+G4double G4hLowEnergyLoss::LOGRTable = 1.1;
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
@@ -772,5 +773,4 @@
     theRangeTable = theRangepbarTable ;
   }
- 
   G4double R2 = RTable*RTable ;
   G4double R1 = RTable+1.;
@@ -1062,5 +1062,6 @@
 {
   G4double LowEdgeRange,A,B,C,discr,KineticEnergy ;
-  G4double Tbin = LowestKineticEnergy/RTable ;
+  G4double Tbin = 0;
+  if (RTable !=0.) Tbin = LowestKineticEnergy/RTable ;
   G4double rangebin = 0.0 ;
   G4int binnumber = -1 ;
