Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundAlpha.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundAlpha.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundAlpha.cc	(revision 1007)
@@ -24,39 +24,27 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundAlpha.cc,v 1.5 2009/02/13 18:57:32 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundAlpha
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada add choice of options  
-// 10.02.2009 J. M. Quesada set default opt1  
-//
+//J.M.Quesada (August 08). New source file
+// 
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
 
 #include "G4PreCompoundAlpha.hh"
 
-G4ReactionProduct * G4PreCompoundAlpha::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct =
-    new G4ReactionProduct(G4Alpha::AlphaDefinition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+  G4ReactionProduct * G4PreCompoundAlpha::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct =
+      new G4ReactionProduct(G4Alpha::AlphaDefinition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}   
-
-G4double G4PreCompoundAlpha::FactorialFactor(const G4double N, const G4double P)
-{
-  return 
+    return theReactionProduct;
+  }   
+
+
+   G4double G4PreCompoundAlpha::FactorialFactor(const G4double N, const G4double P)
+  {
+     return 
       (N-4.0)*(P-3.0)*(
 		       (((N-3.0)*(P-2.0))/2.0) *(
@@ -66,21 +54,23 @@
 						 )
 		       );
-}
-  
-G4double G4PreCompoundAlpha::CoalescenceFactor(const G4double A)
-{
-  return 4096.0/(A*A*A);  
-}    
-
-G4double G4PreCompoundAlpha::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2)*(NumberParticles-3);
-  if(NumberCharged >=2 && (NumberParticles-NumberCharged) >=2 ) {
-    rj = 6.0*static_cast<G4double>(NumberCharged*(NumberCharged-1)*(NumberParticles-NumberCharged)*
-				   (NumberParticles-NumberCharged-1))/static_cast<G4double>(denominator);  
-  }
-  return rj;
-}
+  }
+  
+   G4double G4PreCompoundAlpha::CoalescenceFactor(const G4double A)
+  {
+     return 4096.0/(A*A*A);  
+  }    
+
+
+
+   G4double G4PreCompoundAlpha::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2)*(NumberParticles-3);
+    if(NumberCharged >=2 && (NumberParticles-NumberCharged) >=2 ) rj = 6.0*static_cast<G4double>(NumberCharged*(NumberCharged-1)*(NumberParticles-NumberCharged)*(NumberParticles-NumberCharged-1))/static_cast<G4double>(denominator);  
+ return rj;
+  }
+
+
+
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -90,5 +80,5 @@
 //OPT=3,4 Kalbach's parameterization 
 // 
-G4double G4PreCompoundAlpha::CrossSection(const  G4double K)
+ G4double G4PreCompoundAlpha::CrossSection(const  G4double K)
 {
 
@@ -125,33 +115,33 @@
 //----------------
 //
-G4double G4PreCompoundAlpha::GetAlpha()
-{
-  G4double C = 0.0;
-  G4double aZ = GetZ() + GetRestZ();
-  if (aZ <= 30) 
-    {
-      C = 0.10;
-    } 
-  else if (aZ <= 50) 
-    {
-      C = 0.1 + -((aZ-50.)/20.)*0.02;
-    } 
-  else if (aZ < 70) 
-    {
-      C = 0.08 + -((aZ-70.)/20.)*0.02;
-    }
-  else 
-    {
-      C = 0.06;
-    }
-  return 1.0+C;
-}
+ G4double G4PreCompoundAlpha::GetAlpha()
+  {
+ G4double C = 0.0;
+    G4double aZ = GetZ() + GetRestZ();
+    if (aZ <= 30) 
+      {
+	C = 0.10;
+      } 
+    else if (aZ <= 50) 
+      {
+	C = 0.1 + -((aZ-50.)/20.)*0.02;
+      } 
+    else if (aZ < 70) 
+      {
+	C = 0.08 + -((aZ-70.)/20.)*0.02;
+      }
+    else 
+      {
+	C = 0.06;
+      }
+    return 1.0+C;
+  }
 //
 //--------------------
 // 
-G4double G4PreCompoundAlpha::GetBeta() 
-{
-  return -GetCoulombBarrier();
-}
+  G4double G4PreCompoundAlpha::GetBeta() 
+  {
+      return -GetCoulombBarrier();
+  }
 //
 //********************* OPT=1,2 : Chatterjee's cross section ************************ 
@@ -163,5 +153,5 @@
   G4double Kc=K;
 
-  // JMQ xsec is set constant above limit of validity
+  // JMQ xsec is set constat above limit of validity
   if (K>50) Kc=50;
 
@@ -223,7 +213,5 @@
   G4double      ra=1.20;
         
-  //JMQ 13/02/09 increase of reduced radius to lower the barrier
-  // ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
-  ec = 1.44 * theZ * ResidualZ / (1.7*ResidualAthrd+ra);
+  ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
   ecsq = ec * ec;
   p = p0 + p1/ec + p2/ecsq;
@@ -248,5 +236,5 @@
   
   if (elab <= ec) { //start for E<Ec
-    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;
+    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;    
   }           //end for E<Ec
   else {           //start for E>Ec
@@ -264,2 +252,6 @@
 
 //   ************************** end of cross sections ******************************* 
+
+
+
+
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundDeuteron.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundDeuteron.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundDeuteron.cc	(revision 1007)
@@ -24,56 +24,44 @@
 // ********************************************************************
 //
-// $Id: G4PreCompoundDeuteron.cc,v 1.5 2009/02/13 18:57:32 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundDeuteron
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada add choice of options  
-// 10.02.2009 J. M. Quesada set default opt1  
-//
+//J.M.Quesada (August 08). New source file
+// 
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
 
 #include "G4PreCompoundDeuteron.hh"
 
-G4ReactionProduct * G4PreCompoundDeuteron::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct =
-    new G4ReactionProduct(G4Deuteron::DeuteronDefinition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+
+
+  G4ReactionProduct * G4PreCompoundDeuteron::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct =
+      new G4ReactionProduct(G4Deuteron::DeuteronDefinition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}   
-
- 
-G4double G4PreCompoundDeuteron::FactorialFactor(const G4double N, const G4double P)
-{
-  return (N-1.0)*(N-2.0)*(P-1.0)*P/2.0;
-}
+    return theReactionProduct;
+  }   
+
+ 
+   G4double G4PreCompoundDeuteron::FactorialFactor(const G4double N, const G4double P)
+  {
+      return (N-1.0)*(N-2.0)*(P-1.0)*P/2.0;
+  }
   
-G4double G4PreCompoundDeuteron::CoalescenceFactor(const G4double A)
-{
-  return 16.0/A;
-}    
-
-G4double G4PreCompoundDeuteron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  G4double denominator = NumberParticles*(NumberParticles-1);
-  if(NumberCharged >=1 && (NumberParticles-NumberCharged) >=1) {
-    rj = 2.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged))
-      / static_cast<G4double>(denominator); 
-  }
-  return rj;
-}
+   G4double G4PreCompoundDeuteron::CoalescenceFactor(const G4double A)
+  {
+    return 16.0/A;
+  }    
+
+  G4double G4PreCompoundDeuteron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    G4double denominator = NumberParticles*(NumberParticles-1);
+   if(NumberCharged >=1 && (NumberParticles-NumberCharged) >=1) rj = 2.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged))/static_cast<G4double>(denominator); 
+
+    return rj;
+  }
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -83,5 +71,5 @@
 //OPT=3,4 Kalbach's parameterization 
 // 
-G4double G4PreCompoundDeuteron::CrossSection(const  G4double K)
+ G4double G4PreCompoundDeuteron::CrossSection(const  G4double K)
 {
 
@@ -118,25 +106,25 @@
 //---------
 //
-G4double G4PreCompoundDeuteron::GetAlpha()
-{
-  G4double C = 0.0;
-  G4double aZ = GetZ() + GetRestZ();
-  if (aZ >= 70) 
-    {
-      C = 0.10;
-    } 
-  else 
-    {
-      C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375; 
-    }
-  return 1.0 + C/2.0;
-}
+ G4double G4PreCompoundDeuteron::GetAlpha()
+  {
+G4double C = 0.0;
+    G4double aZ = GetZ() + GetRestZ();
+    if (aZ >= 70) 
+      {
+	C = 0.10;
+      } 
+    else 
+      {
+	C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375; 
+      }
+    return 1.0 + C/2.0;
+  }
 //
 //---------
 //
-G4double G4PreCompoundDeuteron::GetBeta() 
-{
-  return -GetCoulombBarrier();
-}
+  G4double G4PreCompoundDeuteron::GetBeta() 
+  {
+      return -GetCoulombBarrier();
+  }
 //
 //********************* OPT=1,2 : Chatterjee's cross section ************************ 
@@ -152,5 +140,5 @@
 
   G4double landa ,mu ,nu ,p , Ec,q,r,ji,xs;
-  //G4double Eo(0),epsilon1(0),epsilon2(0),discri(0);
+//G4double Eo(0),epsilon1(0),epsilon2(0),discri(0);
 
  
@@ -186,4 +174,7 @@
 }
 
+
+
+
 // *********** OPT=3,4 : Kalbach's cross sections (from PRECO code)*************
 G4double G4PreCompoundDeuteron::GetOpt34(const  G4double K)
@@ -214,7 +205,5 @@
   G4double     ra=0.80;
         
-  //JMQ 13/02/09 increase of reduced radius to lower the barrier
-  // ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
-  ec = 1.44 * theZ * ResidualZ / (1.7*ResidualAthrd+ra);
+  ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
   ecsq = ec * ec;
   p = p0 + p1/ec + p2/ecsq;
@@ -237,7 +226,7 @@
   elab = K * FragmentA / ResidualA;
   sig = 0.;
-
+ 
   if (elab <= ec) { //start for E<Ec
-    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;
+    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;    
   }           //end for E<Ec
   else {           //start for E>Ec
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundEmission.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundEmission.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundEmission.cc	(revision 1007)
@@ -25,18 +25,7 @@
 //
 //
-// $Id: G4PreCompoundEmission.cc,v 1.19 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundEmission
-//
-// Author:         V.Lara
-//
-// Modified:  
-//
+// Hadronic Process: Nuclear Preequilibrium
+// by V. Lara 
+
 
 #include "G4PreCompoundEmission.hh"
@@ -53,10 +42,10 @@
 
 
-G4bool G4PreCompoundEmission::operator==(const G4PreCompoundEmission &) const
+  G4bool G4PreCompoundEmission::operator==(const G4PreCompoundEmission &) const
 {
   return false;
 }
 
-G4bool G4PreCompoundEmission::operator!=(const G4PreCompoundEmission &) const
+    G4bool G4PreCompoundEmission::operator!=(const G4PreCompoundEmission &) const
 {
   return true;
@@ -106,4 +95,6 @@
   return;
 }
+
+
 
 G4ReactionProduct * G4PreCompoundEmission::PerformEmission(G4Fragment & aFragment)
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragment.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragment.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragment.cc	(revision 1007)
@@ -24,9 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PreCompoundFragment.cc,v 1.8 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// J. M. Quesada (August 2008).  
-// Based  on previous work by V. Lara
+//J. M. Quesada (August 2008).  
+//Based  on previous work by V. Lara
 // JMQ (06 September 2008) Also external choice has been added for:
 //                      - superimposed Coulomb barrier (if useSICB=true) 
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragmentVector.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragmentVector.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundFragmentVector.cc	(revision 1007)
@@ -24,6 +24,7 @@
 // ********************************************************************
 //
-// $Id: G4PreCompoundFragmentVector.cc,v 1.11 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+//
+// $Id: G4PreCompoundFragmentVector.cc,v 1.7 2008/05/08 10:36:03 quesada Exp $
+// GEANT4 tag $Name: geant4-09-02 $
 //
 // Hadronic Process: Nuclear Preequilibrium
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundHe3.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundHe3.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundHe3.cc	(revision 1007)
@@ -24,39 +24,28 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundHe3.cc,v 1.5 2009/02/13 18:57:32 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundHe3
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada add choice of options  
-// 10.02.2009 J. M. Quesada set default opt1  
-//
+//J.M.Quesada (August 08). New source file
+// 
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
  
 #include "G4PreCompoundHe3.hh"
 
-G4ReactionProduct * G4PreCompoundHe3::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct =
-    new G4ReactionProduct(G4He3::He3Definition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+
+  G4ReactionProduct * G4PreCompoundHe3::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct =
+      new G4ReactionProduct(G4He3::He3Definition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}   
-
-G4double G4PreCompoundHe3::FactorialFactor(const G4double N, const G4double P)
-{
-  return 
+    return theReactionProduct;
+  }   
+
+
+   G4double G4PreCompoundHe3::FactorialFactor(const G4double N, const G4double P)
+  {
+     return 
       (N-3.0)*(P-2.0)*(
 		       (((N-2.0)*(P-1.0))/2.0) *(
@@ -64,21 +53,24 @@
 						 )
 		       );
-}
-  
-G4double G4PreCompoundHe3::CoalescenceFactor(const G4double A)
-{
-  return 243.0/(A*A);
-}    
-
-G4double G4PreCompoundHe3::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2);
-  if(NumberCharged >=2 && (NumberParticles-NumberCharged) >= 1) {
-    rj = 3.0*static_cast<G4double>(NumberCharged*(NumberCharged-1)*(NumberParticles-NumberCharged))
-      / static_cast<G4double>(denominator);  
-  }
-  return rj;
-}
+  }
+  
+  G4double G4PreCompoundHe3::CoalescenceFactor(const G4double A)
+  {
+         return 243.0/(A*A);
+  }    
+
+
+
+  G4double G4PreCompoundHe3::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2);
+    if(NumberCharged >=2 && (NumberParticles-NumberCharged) >= 1) rj = 3.0*static_cast<G4double>(NumberCharged*(NumberCharged-1)*(NumberParticles-NumberCharged))/static_cast<G4double>(denominator);  
+    
+    return rj;
+  }
+
+
+
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -88,5 +80,5 @@
 //OPT=3,4 Kalbach's parameterization 
 // 
-G4double G4PreCompoundHe3::CrossSection(const  G4double K)
+ G4double G4PreCompoundHe3::CrossSection(const  G4double K)
 {
   ResidualA=GetRestA();
@@ -122,33 +114,33 @@
 //----------------
 //
-G4double G4PreCompoundHe3::GetAlpha()
-{
-  G4double C = 0.0;
-  G4double aZ = GetZ() + GetRestZ();
-  if (aZ <= 30) 
-    {
-      C = 0.10;
-    }
-  else if (aZ <= 50) 
-    {
-      C = 0.1 + -((aZ-50.)/20.)*0.02;
-    } 
-  else if (aZ < 70) 
-    {
-      C = 0.08 + -((aZ-70.)/20.)*0.02;
-    }
-  else 
-    {
-      C = 0.06;
-    }
-  return 1.0 + C*(4.0/3.0);
-}
+  G4double G4PreCompoundHe3::GetAlpha()
+  {
+    G4double C = 0.0;
+    G4double aZ = GetZ() + GetRestZ();
+    if (aZ <= 30) 
+      {
+	C = 0.10;
+      }
+    else if (aZ <= 50) 
+      {
+	C = 0.1 + -((aZ-50.)/20.)*0.02;
+      } 
+    else if (aZ < 70) 
+      {
+	C = 0.08 + -((aZ-70.)/20.)*0.02;
+      }
+    else 
+      {
+	C = 0.06;
+      }
+    return 1.0 + C*(4.0/3.0);
+  }
 //
 //--------------------
 // 
-G4double G4PreCompoundHe3::GetBeta() 
-{
-  return -GetCoulombBarrier();
-}
+   G4double G4PreCompoundHe3::GetBeta() 
+  {
+      return -GetCoulombBarrier();
+  }
 //
 //********************* OPT=1,2 : Chatterjee's cross section ************************ 
@@ -222,7 +214,5 @@
   G4double      ra=0.80;
         
-  //JMQ 13/02/09 increase of reduced radius to lower the barrier
-  // ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
-  ec = 1.44 * theZ * ResidualZ / (1.7*ResidualAthrd+ra);
+  ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
   ecsq = ec * ec;
   p = p0 + p1/ec + p2/ecsq;
@@ -247,5 +237,5 @@
   
   if (elab <= ec) { //start for E<Ec
-    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;
+    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;    
   }           //end for E<Ec
   else {           //start for E>Ec
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundIon.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundIon.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundIon.cc	(revision 1007)
@@ -24,18 +24,6 @@
 // ********************************************************************
 //
-// $Id: G4PreCompoundIon.cc,v 1.16 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundIon
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 10.02.2009 J. M. Quesada fixed bug in density level of light fragments  
+//J. M. Quesada (August 2008).  
+//Based  on previous work by V. Lara
 //
 
@@ -43,10 +31,10 @@
 #include "G4PreCompoundParameters.hh"
 
-G4bool G4PreCompoundIon::IsItPossible(const G4Fragment& aFragment) 
-{
-  G4int pplus = aFragment.GetNumberOfCharged();   
-  G4int pneut = aFragment.GetNumberOfParticles()-pplus;
-  return (pneut >= (GetA()-GetZ()) && pplus >= GetZ());
-}
+ G4bool G4PreCompoundIon::IsItPossible(const G4Fragment& aFragment) 
+  {
+    G4int pplus = aFragment.GetNumberOfCharged();   
+    G4int pneut = aFragment.GetNumberOfParticles()-pplus;
+    return (pneut >= (GetA()-GetZ()) && pplus >= GetZ());
+  }
 
 G4double G4PreCompoundIon::
@@ -69,9 +57,6 @@
     G4PreCompoundParameters::GetAddress()->GetLevelDensity();
 
-  //JMQ 06/02/209  This is  THE BUG that was killing cluster emission
-  // G4double gj = (6.0/pi2)*GetA() *
-  //   G4PreCompoundParameters::GetAddress()->GetLevelDensity();
-
-  G4double gj = g1;
+  G4double gj = (6.0/pi2)*GetA() *
+    G4PreCompoundParameters::GetAddress()->GetLevelDensity();
 
   G4double A0 = ((P*P+H*H+P-H)/4.0 - H/2.0)/g0; 
@@ -89,10 +74,7 @@
   G4double Ej = std::max(0.0,eKin + GetBindingEnergy() -Aj); 
 
-  // JMQ 10/02/09 reshaping of the formula (unnecessary std::pow elimitated)
-  G4double pA = (3.0/4.0) * std::sqrt(std::max(0.0, 2.0/(GetReducedMass()*
-		(eKin+GetBindingEnergy()))))/(pi * r0 * r0 *r0* GetRestA())* 
-                eKin*CrossSection(eKin)*millibarn* 
-                CoalescenceFactor(aFragment.GetA()) * FactorialFactor(N,P)* 
-                GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged());
+
+ G4double pA = 1.e-25*(3.0/4.0) * std::sqrt(std::max(0.0, 2.0/(GetReducedMass()*
+(eKin+GetBindingEnergy()))))/(pi * r0 * r0 * std::pow(GetRestA(),2.0/3.0) )* eKin*CrossSection(eKin) /(r0*std::pow(GetRestA(),1.0/3.0)) * CoalescenceFactor(aFragment.GetA()) * FactorialFactor(N,P)* GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())  ;
 
   G4double pB = std::pow((g1*E1)/(g0*E0),N-GetA()-1.0)*(g1/g0);
@@ -100,6 +82,13 @@
   G4double pC = std::pow((gj*Ej)/(g0*E0),GetA()-1.0)*(gj/g0)/E0; 
 
+
   G4double Probability = pA * pB * pC;
 
   return Probability;
 }
+
+
+
+
+
+
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundModel.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundModel.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundModel.cc	(revision 1007)
@@ -26,5 +26,5 @@
 //
 // $Id: G4PreCompoundModel.cc,v 1.17 2008/12/09 14:09:59 ahoward Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+// GEANT4 tag $Name: geant4-09-02 $
 //
 // by V. Lara
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNeutron.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNeutron.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNeutron.cc	(revision 1007)
@@ -24,43 +24,31 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundNeutron.cc,v 1.4 2009/02/11 18:06:00 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundNeutron
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada add choice of options  
-// 10.02.2009 J. M. Quesada set default opt3
+//J. M. Quesada (August 2008). New source file
+//
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
 // 
-
 #include "G4PreCompoundNeutron.hh"
 
-G4ReactionProduct * G4PreCompoundNeutron::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct = 
-    new G4ReactionProduct(G4Neutron::NeutronDefinition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+
+  G4ReactionProduct * G4PreCompoundNeutron::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct = 
+      new G4ReactionProduct(G4Neutron::NeutronDefinition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}
-
-G4double G4PreCompoundNeutron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  if(NumberParticles > 0) rj = static_cast<G4double>(NumberParticles - NumberCharged)/
-			    static_cast<G4double>(NumberParticles);
-  return rj;
-}
+    return theReactionProduct;
+  }
+
+
+   G4double G4PreCompoundNeutron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    if(NumberParticles > 0) rj = static_cast<G4double>(NumberParticles - NumberCharged)/static_cast<G4double>(NumberParticles);
+    return rj;
+  }
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -70,5 +58,5 @@
 //OPT=3,4 Kalbach's parameterization 
 // 
-G4double G4PreCompoundNeutron::CrossSection(const  G4double K)
+ G4double G4PreCompoundNeutron::CrossSection(const  G4double K)
 {
   ResidualA=GetRestA();
@@ -104,18 +92,21 @@
 //-------
 //
-G4double G4PreCompoundNeutron::GetAlpha()
-{
-  //   return 0.76+2.2/std::pow(GetRestA(),1.0/3.0);
+  G4double G4PreCompoundNeutron::GetAlpha()
+  {
+ //   return 0.76+2.2/std::pow(GetRestA(),1.0/3.0);
   return 0.76+2.2/ResidualAthrd;
-}
+  }
 //
 //------------
 //
-G4double G4PreCompoundNeutron::GetBeta() 
-{
-  //   return (2.12/std::pow(GetRestA(),2.0/3.0)-0.05)*MeV/GetAlpha();
-  return (2.12/(ResidualAthrd*ResidualAthrd)-0.05)*MeV/GetAlpha();
-}
-//
+  G4double G4PreCompoundNeutron::GetBeta() 
+  {
+ //   return (2.12/std::pow(GetRestA(),2.0/3.0)-0.05)*MeV/GetAlpha();
+ return (2.12/(ResidualAthrd*ResidualAthrd)-0.05)*MeV/GetAlpha();
+  }
+//
+
+  
+
 
 //********************* OPT=1,2 : Chatterjee's cross section ************************ 
@@ -165,5 +156,5 @@
   G4double b,ecut,cut,ecut2,geom,elab;
 
-  //safety initialization
+//safety initialization
   landa0=0;
   landa1=0;
@@ -181,5 +172,7 @@
   spill= 1.0e+18; 
 
-  // PRECO xs for neutrons is choosen
+
+
+// PRECO xs for neutrons is choosen
 
   p0 = -312.;
@@ -208,7 +201,9 @@
   xnulam = 1.;
   etest = 32.;
-  //          ** etest is the energy above which the rxn cross section is
-  //          ** compared with the geometrical limit and the max taken.
-  //          ** xnulam here is a dummy value to be used later.
+//          ** etest is the energy above which the rxn cross section is
+//          ** compared with the geometrical limit and the max taken.
+//          ** xnulam here is a dummy value to be used later.
+
+
 
   a = -2.*p*ec + landa - nu/ecsq;
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc	(revision 1007)
@@ -24,19 +24,6 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundNucleon.cc,v 1.13 2009/02/11 18:06:00 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundNucleon
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 10.02.2009 J. M. Quesada cleanup
+//J. M. Quesada (August 2008).  
+//Based  on previous work by V. Lara
 //
 
@@ -45,9 +32,10 @@
 
 G4bool G4PreCompoundNucleon::IsItPossible(const G4Fragment& aFragment) 
-{
-  G4int pplus = aFragment.GetNumberOfCharged();   
-  G4int pneut = aFragment.GetNumberOfParticles()-pplus;
-  return (pneut >= (GetA()-GetZ()) && pplus >= GetZ());
-}
+  {
+    G4int pplus = aFragment.GetNumberOfCharged();   
+    G4int pneut = aFragment.GetNumberOfParticles()-pplus;
+    return (pneut >= (GetA()-GetZ()) && pplus >= GetZ());
+  }
+
 
 G4double G4PreCompoundNucleon::
@@ -57,4 +45,5 @@
   if ( !IsItPossible(aFragment) ) return 0.0;
 
+
   G4double U = aFragment.GetExcitationEnergy();
   G4double P = aFragment.GetNumberOfParticles();
@@ -62,4 +51,6 @@
   G4double N = P + H;
   
+
+
   G4double g0 = (6.0/pi2)*aFragment.GetA() * 
     G4PreCompoundParameters::GetAddress()->GetLevelDensity();
@@ -67,4 +58,6 @@
   G4double g1 = (6.0/pi2)*GetRestA() * 
     G4PreCompoundParameters::GetAddress()->GetLevelDensity();
+
+
 
   G4double A0 = ((P*P+H*H+P-H)/4.0 - H/2.0)/g0;
@@ -78,18 +71,16 @@
 
 
-  G4double Probability = 1.0/pi2*2.0/(hbarc*hbarc*hbarc) * GetReducedMass()
-    * GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())  
-    * eKin*CrossSection(eKin)*millibarn * P*(N-1.0) * std::pow(g1*E1/(g0*E0),N-2.0)/E0 * g1/(g0*g0);
+ G4double Probability =1.e-25* 1.0/pi2*2.0/(hbarc*hbarc*hbarc) * GetReducedMass()
+* GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())  * eKin*CrossSection(eKin) * P*(N-1.0) * std::pow(g1*E1/(g0*E0),N-2.0)/E0 * g1/(g0*g0);
 
-  if (GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())<0.0 
-      || CrossSection(eKin) <0) {  
-    std::ostringstream errOs;
-    G4cout<<"WARNING:  NEGATIVE VALUES "<<G4endl;     
-    errOs << "Rj=" << GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())
-	  <<G4endl;
-    errOs <<"  xsec("<<eKin<<" MeV) ="<<CrossSection(eKin)<<G4endl;
-    errOs <<"  A="<<GetA()<<"  Z="<<GetZ()<<G4endl;
-    throw G4HadronicException(__FILE__, __LINE__, errOs.str());
-  }
+
+if (GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())<0.0 || CrossSection(eKin) <0)
+{  std::ostringstream errOs;
+   G4cout<<"WARNING:  NEGATIVE VALUES "<<G4endl;     
+        errOs << "Rj=" << GetRj(aFragment.GetNumberOfParticles(), aFragment.GetNumberOfCharged())<<G4endl;
+        errOs <<"  xsec("<<eKin<<" MeV) ="<<CrossSection(eKin)<<G4endl;
+        errOs <<"  A="<<GetA()<<"  Z="<<GetZ()<<G4endl;
+        throw G4HadronicException(__FILE__, __LINE__, errOs.str());
+              }
 
   return Probability;
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundParameters.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundParameters.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundParameters.cc	(revision 1007)
@@ -26,5 +26,5 @@
 //
 // $Id: G4PreCompoundParameters.cc,v 1.3 2006/06/29 20:59:29 gunter Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+// GEANT4 tag $Name: geant4-09-02 $
 //
 // by V. Lara
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundProton.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundProton.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundProton.cc	(revision 1007)
@@ -24,43 +24,35 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundProton.cc,v 1.4 2009/02/11 18:06:00 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundProton
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada added external choice of inverse cross section option
-// 21.08.2008 J. M. Quesada added external choice for superimposed Coulomb barrier 
-//                          (if useSICB=true) 
-//
+//J.M.Quesada (August 08). New source file
+// 
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
+// 
+// JMQ (06 September 2008) Also external choice has been added for:
+//                      - superimposed Coulomb barrier (if useSICB=true) 
 
 #include "G4PreCompoundProton.hh"
 
-G4ReactionProduct * G4PreCompoundProton::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct = 
-    new G4ReactionProduct(G4Proton::ProtonDefinition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+  G4ReactionProduct * G4PreCompoundProton::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct = 
+      new G4ReactionProduct(G4Proton::ProtonDefinition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}
-
-G4double G4PreCompoundProton::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  if(NumberParticles > 0) rj = static_cast<G4double>(NumberCharged)/static_cast<G4double>(NumberParticles);
-  return rj;
-}
+    return theReactionProduct;
+  }
+
+
+   G4double G4PreCompoundProton::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    if(NumberParticles > 0) rj = static_cast<G4double>(NumberCharged)/static_cast<G4double>(NumberParticles);
+    return rj;
+  }
+
+
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -71,5 +63,5 @@
 //OPT=3 Kalbach's parameterization 
 // 
-G4double G4PreCompoundProton::CrossSection(const  G4double K)
+ G4double G4PreCompoundProton::CrossSection(const  G4double K)
 {
   //G4cout<<" In G4PreCompoundProton OPTxs="<<OPTxs<<G4endl;
@@ -83,4 +75,6 @@
   FragmentA=GetA()+GetRestA();
   FragmentAthrd=std::pow(FragmentA,0.33333);
+
+
 
   if (OPTxs==0) return GetOpt0(K);
@@ -100,35 +94,37 @@
 G4double G4PreCompoundProton::GetOpt0(const  G4double K)
 {
-  const G4double r0 = G4PreCompoundParameters::GetAddress()->Getr0();
+const G4double r0 = G4PreCompoundParameters::GetAddress()->Getr0();
   // cross section is now given in mb (r0 is in mm) for the sake of consistency
   //with the rest of the options
-  return 1.e+25*pi*(r0*ResidualAthrd)*(r0*ResidualAthrd)*GetAlpha()*(1.+GetBeta()/K);
+ return 1.e+25*pi*(r0*ResidualAthrd)*(r0*ResidualAthrd)*GetAlpha()*(1.+GetBeta()/K);
 }
 //
 //------------
 //
-G4double G4PreCompoundProton::GetAlpha()
-{
-  G4double aZ = static_cast<G4double>(GetRestZ());
-  G4double C = 0.0;
-  if (aZ >= 70) 
-    {
-      C = 0.10;
-    } 
-  else 
-    {
-      C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375;
-    }
-  return 1.0 + C;
-}
+  G4double G4PreCompoundProton::GetAlpha()
+  {
+ G4double aZ = static_cast<G4double>(GetRestZ());
+    G4double C = 0.0;
+    if (aZ >= 70) 
+      {
+	C = 0.10;
+      } 
+    else 
+      {
+	C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375;
+      }
+    return 1.0 + C;
+  }
 //
 //-------------------
 //  
-G4double G4PreCompoundProton::GetBeta() 
-{
+  G4double G4PreCompoundProton::GetBeta() 
+  {
   return -GetCoulombBarrier();
-}
-//
-  
+  }
+//
+  
+
+
 //********************* OPT=1 : Chatterjee's cross section ************************ 
 //(fitting to cross section from Bechetti & Greenles OM potential)
@@ -136,40 +132,42 @@
 G4double G4PreCompoundProton::GetOpt1(const  G4double K)
 {
-  G4double Kc=K; 
-
-  // JMQ  xsec is set constat above limit of validity
-  if (K>50)  Kc=50;
-
-  G4double landa, landa0, landa1, mu, mu0, mu1,nu, nu0, nu1, nu2,xs;
-  G4double p, p0, p1, p2,Ec,delta,q,r,ji;
-  
-  p0 = 15.72;
-  p1 = 9.65;
-  p2 = -449.0;
-  landa0 = 0.00437;
-  landa1 = -16.58;
-  mu0 = 244.7;
-  mu1 = 0.503;
-  nu0 = 273.1;
-  nu1 = -182.4;
-  nu2 = -1.872;  
-  delta=0.;  
-
-  Ec = 1.44*theZ*ResidualZ/(1.5*ResidualAthrd+delta);
-  p = p0 + p1/Ec + p2/(Ec*Ec);
-  landa = landa0*ResidualA + landa1;
-  mu = mu0*std::pow(ResidualA,mu1);
-  nu = std::pow(ResidualA,mu1)*(nu0 + nu1*Ec + nu2*(Ec*Ec));
-  q = landa - nu/(Ec*Ec) - 2*p*Ec;
-  r = mu + 2*nu/Ec + p*(Ec*Ec);
-
-  ji=std::max(Kc,Ec);
-  if(Kc < Ec) { xs = p*Kc*Kc + q*Kc + r;}
-  else {xs = p*(Kc - ji)*(Kc - ji) + landa*Kc + mu + nu*(2 - Kc/ji)/ji ;}
-  if (xs <0.0) {xs=0.0;}
-
-  return xs; 
+ G4double Kc=K; 
+
+// JMQ  xsec is set constat above limit of validity
+ if (K>50)  Kc=50;
+
+ G4double landa, landa0, landa1, mu, mu0, mu1,nu, nu0, nu1, nu2,xs;
+ G4double p, p0, p1, p2,Ec,delta,q,r,ji;
+  
+ p0 = 15.72;
+ p1 = 9.65;
+ p2 = -449.0;
+ landa0 = 0.00437;
+ landa1 = -16.58;
+ mu0 = 244.7;
+ mu1 = 0.503;
+ nu0 = 273.1;
+ nu1 = -182.4;
+ nu2 = -1.872;  
+ delta=0.;  
+
+ Ec = 1.44*theZ*ResidualZ/(1.5*ResidualAthrd+delta);
+ p = p0 + p1/Ec + p2/(Ec*Ec);
+ landa = landa0*ResidualA + landa1;
+ mu = mu0*std::pow(ResidualA,mu1);
+ nu = std::pow(ResidualA,mu1)*(nu0 + nu1*Ec + nu2*(Ec*Ec));
+ q = landa - nu/(Ec*Ec) - 2*p*Ec;
+ r = mu + 2*nu/Ec + p*(Ec*Ec);
+
+ ji=std::max(Kc,Ec);
+ if(Kc < Ec) { xs = p*Kc*Kc + q*Kc + r;}
+ else {xs = p*(Kc - ji)*(Kc - ji) + landa*Kc + mu + nu*(2 - Kc/ji)/ji ;}
+ if (xs <0.0) {xs=0.0;}
+
+ return xs; 
 
 }
+
+
 
 //************* OPT=2 : Welisch's proton reaction cross section ************************ 
@@ -180,6 +178,6 @@
   G4double rnpro,rnneu,eekin,ekin,ff1,ff2,ff3,r0,fac,fac1,fac2,b0,xine_th(0);
  
-  //This is redundant when the Coulomb  barrier is overimposed to all cross sections 
-  //It should be kept when Coulomb barrier only imposed at OPTxs=2
+//This is redundant when the Coulomb  barrier is overimposed to all cross sections 
+//It should be kept when Coulomb barrier only imposed at OPTxs=2
 
   if(!useSICB && K<=theCoulombBarrier) return xine_th=0.0;
@@ -318,6 +316,7 @@
   }   //end for E>Ec
 
-  return sig;
-}
+  return sig;}
+
+
 
 //   ************************** end of cross sections ******************************* 
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTransitions.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTransitions.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTransitions.cc	(revision 1007)
@@ -24,22 +24,10 @@
 // ********************************************************************
 //
-// $Id: G4PreCompoundTransitions.cc,v 1.20 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
 //
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundIon
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 16.02.2008 J. M. Quesada fixed bugs 
-// 06.09.2008 J. M. Quesada added external choices for:
+//J. M. Quesada (Feb. 08). Base on previous work by V. Lara.
+// New transition probabilities. Several bugs fixed.
+// JMQ (06 September 2008) Also external choices have been added for:
 //                      - "never go back"  hipothesis (useNGB=true) 
 //                      -  CEM transition probabilities (useCEMtr=true)
-
 #include "G4PreCompoundTransitions.hh"
 #include "G4HadronicException.hh"
@@ -196,6 +184,5 @@
     //  G4cout<<"N="<<N<<"  P="<<P<<"  H="<<H<<G4endl;
     //  G4cout<<"l+ ="<<TransitionProb1<<"  l- ="<< TransitionProb2<<"  l0 ="<< TransitionProb3<<G4endl; 
-    return TransitionProb1 + TransitionProb2 + TransitionProb3;
-  }
+    return TransitionProb1 + TransitionProb2 + TransitionProb3;}
   
   else  {
@@ -228,4 +215,6 @@
     return TransitionProb1 + TransitionProb2 + TransitionProb3;
   }
+
+
 }
 
@@ -248,21 +237,17 @@
     }
 
-  // AH/JMQ: Randomly decrease the number of charges if deltaN is -2 and in proportion 
-  // to the number charges w.r.t. number of particles,  PROVIDED that there are charged particles
-  if(deltaN < 0 && G4UniformRand() <= 
-     static_cast<G4double>(result.GetNumberOfCharged())/static_cast<G4double>(result.GetNumberOfParticles()) 
-     && (result.GetNumberOfCharged() >= 1)) {
+  // AH/JMQ: Randomly decrease the number of charges if deltaN is -2 and in proportion to the number charges w.r.t. number of particles,  PROVIDED that there are charged particles
+  if(deltaN < 0 && G4UniformRand() <= static_cast<G4double>(result.GetNumberOfCharged())/static_cast<G4double>(result.GetNumberOfParticles()) && (result.GetNumberOfCharged() >= 1))
     result.SetNumberOfCharged(result.GetNumberOfCharged()+deltaN/2); // deltaN is negative!
-  }
-
-  // JMQ the following lines have to be before SetNumberOfCharged, otherwise the check on 
-  // number of charged vs. number of particles fails
+
+
+ //JMQ the following lines have to be before SetNumberOfCharged, otherwise the check on number of charged vs. number of particles fails
   result.SetNumberOfParticles(result.GetNumberOfParticles()+deltaN/2);
   result.SetNumberOfHoles(result.GetNumberOfHoles()+deltaN/2); 
 
-  // With weight Z/A, number of charged particles is increased with +1
+ // With weight Z/A, number of charged particles is increased with +1
   if ( ( deltaN > 0 ) &&
       (G4UniformRand() <= static_cast<G4double>(result.GetZ()-result.GetNumberOfCharged())/
-       std::max(static_cast<G4double>(result.GetA()-Nexcitons),1.)))
+		  std::max(static_cast<G4double>(result.GetA()-Nexcitons),1.)))
     {
       result.SetNumberOfCharged(result.GetNumberOfCharged()+deltaN/2);
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTriton.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTriton.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundTriton.cc	(revision 1007)
@@ -25,39 +25,27 @@
 //
 //
-// $Id: G4PreCompoundTriton.cc,v 1.5 2009/02/13 18:57:32 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
-//
-// -------------------------------------------------------------------
-//
-// GEANT4 Class file
-//
-//
-// File name:     G4PreCompoundTriton
-//
-// Author:         V.Lara
-//
-// Modified:  
-// 21.08.2008 J. M. Quesada add choice of options  
-// 10.02.2009 J. M. Quesada set default opt1  
-//
+//J.M.Quesada (August 08). New source file
+//
+// Modif (21 August 2008) by J. M. Quesada for external choice of inverse 
+// cross section option
  
 #include "G4PreCompoundTriton.hh"
 
 
-G4ReactionProduct * G4PreCompoundTriton::GetReactionProduct() const
-{
-  G4ReactionProduct * theReactionProduct =
-    new G4ReactionProduct(G4Triton::TritonDefinition());
-  theReactionProduct->SetMomentum(GetMomentum().vect());
-  theReactionProduct->SetTotalEnergy(GetMomentum().e());
+  G4ReactionProduct * G4PreCompoundTriton::GetReactionProduct() const
+  {
+    G4ReactionProduct * theReactionProduct =
+      new G4ReactionProduct(G4Triton::TritonDefinition());
+    theReactionProduct->SetMomentum(GetMomentum().vect());
+    theReactionProduct->SetTotalEnergy(GetMomentum().e());
 #ifdef PRECOMPOUND_TEST
-  theReactionProduct->SetCreatorModel("G4PrecompoundModel");
+    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
 #endif
-  return theReactionProduct;
-}   
-
-G4double G4PreCompoundTriton::FactorialFactor(const G4double N, const G4double P)
-{
-  return 
+    return theReactionProduct;
+  }   
+
+   G4double G4PreCompoundTriton::FactorialFactor(const G4double N, const G4double P)
+  {
+      return 
       (N-3.0)*(P-2.0)*(
 		       (((N-2.0)*(P-1.0))/2.0) *(
@@ -65,21 +53,22 @@
 						 )
 		       );
-}
+  }
   
-G4double G4PreCompoundTriton::CoalescenceFactor(const G4double A)
-{
-  return 243.0/(A*A);
-}    
-
-G4double G4PreCompoundTriton::GetRj(const G4int NumberParticles, const G4int NumberCharged)
-{
-  G4double rj = 0.0;
-  G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2);
-  if(NumberCharged >= 1 && (NumberParticles-NumberCharged) >= 2) {
-    rj = 3.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged)*(NumberParticles-NumberCharged-1))
-      /static_cast<G4double>(denominator); 
-  }
-  return rj;
-}
+   G4double G4PreCompoundTriton::CoalescenceFactor(const G4double A)
+  {
+     return 243.0/(A*A);
+  }    
+
+
+   G4double G4PreCompoundTriton::GetRj(const G4int NumberParticles, const G4int NumberCharged)
+  {
+    G4double rj = 0.0;
+    G4double denominator = NumberParticles*(NumberParticles-1)*(NumberParticles-2);
+    if(NumberCharged >= 1 && (NumberParticles-NumberCharged) >= 2) rj = 3.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged)*(NumberParticles-NumberCharged-1))/static_cast<G4double>(denominator); 
+    return rj;
+  }
+
+
+
 
 ////////////////////////////////////////////////////////////////////////////////////
@@ -89,5 +78,5 @@
 //OPT=3,4 Kalbach's parameterization 
 // 
-G4double G4PreCompoundTriton::CrossSection(const  G4double K)
+ G4double G4PreCompoundTriton::CrossSection(const  G4double K)
 {
   ResidualA=GetRestA();
@@ -122,26 +111,26 @@
 //---------
 //
-G4double G4PreCompoundTriton::GetAlpha()
-{
-  G4double C = 0.0;
-  G4double aZ = GetZ() + GetRestZ();
-  if (aZ >= 70) 
-    {
-      C = 0.10;
-    } 
-  else 
-    {
-      C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375; 
-    }
- 
-  return 1.0 + C/3.0;
-}
+  G4double G4PreCompoundTriton::GetAlpha()
+  {
+    G4double C = 0.0;
+    G4double aZ = GetZ() + GetRestZ();
+    if (aZ >= 70) 
+      {
+	C = 0.10;
+      } 
+    else 
+      {
+	C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375; 
+      }
+ 
+    return 1.0 + C/3.0;
+  }
 //
 //-------------
 //
-G4double G4PreCompoundTriton::GetBeta() 
-{
-  return -GetCoulombBarrier();
-}
+   G4double G4PreCompoundTriton::GetBeta() 
+  {
+      return -GetCoulombBarrier();
+  }
 //
 //********************* OPT=1,2 : Chatterjee's cross section ************************ 
@@ -194,64 +183,66 @@
 
   G4double landa, mu, nu, p , signor(1.),sig;
-  G4double ec,ecsq,xnulam,etest(0.),a; 
-  G4double b,ecut,cut,ecut2,geom,elab;
-
-
-  G4double     flow = 1.e-18;
-  G4double     spill= 1.e+18;
-
-
-  G4double     p0 = -21.45;
-  G4double     p1 = 484.7;
-  G4double     p2 = -1608.;
-  G4double     landa0 = 0.0186;
-  G4double     landa1 = -8.90;
-  G4double     mu0 = 686.3;
-  G4double     mu1 = 0.325;
-  G4double     nu0 = 368.9;
-  G4double     nu1 = -522.2;
-  G4double     nu2 = -4.998;  
+G4double ec,ecsq,xnulam,etest(0.),a; 
+G4double b,ecut,cut,ecut2,geom,elab;
+
+
+ G4double     flow = 1.e-18;
+ G4double     spill= 1.e+18;
+
+
+ G4double     p0 = -21.45;
+ G4double     p1 = 484.7;
+ G4double     p2 = -1608.;
+ G4double     landa0 = 0.0186;
+ G4double     landa1 = -8.90;
+ G4double     mu0 = 686.3;
+ G4double     mu1 = 0.325;
+ G4double     nu0 = 368.9;
+ G4double     nu1 = -522.2;
+ G4double     nu2 = -4.998;  
   
-  G4double      ra=0.80;
+ G4double      ra=0.80;
         
-  //JMQ 13/02/09 increase of reduced radius to lower the barrier
-  // ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
-  ec = 1.44 * theZ * ResidualZ / (1.7*ResidualAthrd+ra);
-  ecsq = ec * ec;
-  p = p0 + p1/ec + p2/ecsq;
-  landa = landa0*ResidualA + landa1;
-  a = std::pow(ResidualA,mu1);
-  mu = mu0 * a;
-  nu = a* (nu0+nu1*ec+nu2*ecsq);  
-  xnulam = nu / landa;
-  if (xnulam > spill) xnulam=0.;
-  if (xnulam >= flow) etest = 1.2 *std::sqrt(xnulam);
- 
-  a = -2.*p*ec + landa - nu/ecsq;
-  b = p*ecsq + mu + 2.*nu/ec;
-  ecut = 0.;
-  cut = a*a - 4.*p*b;
-  if (cut > 0.) ecut = std::sqrt(cut);
-  ecut = (ecut-a) / (p+p);
-  ecut2 = ecut;
-  if (cut < 0.) ecut2 = ecut - 2.;
-  elab = K * FragmentA / ResidualA;
-  sig = 0.;
- 
-  if (elab <= ec) { //start for E<Ec
-    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;
-  }           //end for E<Ec
-  else {           //start for E>Ec
-    sig = (landa*elab+mu+nu/elab) * signor;
-    geom = 0.;
-    if (xnulam < flow || elab < etest) return sig;
-    geom = std::sqrt(theA*K);
-    geom = 1.23*ResidualAthrd + ra + 4.573/geom;
-    geom = 31.416 * geom * geom;
-    sig = std::max(geom,sig);
-  }           //end for E>Ec
-  return sig;
+ ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
+ ecsq = ec * ec;
+ p = p0 + p1/ec + p2/ecsq;
+ landa = landa0*ResidualA + landa1;
+ a = std::pow(ResidualA,mu1);
+ mu = mu0 * a;
+ nu = a* (nu0+nu1*ec+nu2*ecsq);  
+ xnulam = nu / landa;
+ if (xnulam > spill) xnulam=0.;
+ if (xnulam >= flow) etest = 1.2 *std::sqrt(xnulam);
+ 
+ a = -2.*p*ec + landa - nu/ecsq;
+ b = p*ecsq + mu + 2.*nu/ec;
+ ecut = 0.;
+ cut = a*a - 4.*p*b;
+ if (cut > 0.) ecut = std::sqrt(cut);
+ ecut = (ecut-a) / (p+p);
+ ecut2 = ecut;
+ if (cut < 0.) ecut2 = ecut - 2.;
+ elab = K * FragmentA / ResidualA;
+ sig = 0.;
+ 
+ if (elab <= ec) { //start for E<Ec
+   if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;    
+ }           //end for E<Ec
+ else {           //start for E>Ec
+   sig = (landa*elab+mu+nu/elab) * signor;
+   geom = 0.;
+   if (xnulam < flow || elab < etest) return sig;
+   geom = std::sqrt(theA*K);
+   geom = 1.23*ResidualAthrd + ra + 4.573/geom;
+   geom = 31.416 * geom * geom;
+   sig = std::max(geom,sig);
+ }           //end for E>Ec
+ return sig;
 
 }
 
 //   ************************** end of cross sections ******************************* 
+
+
+
+
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundFragment.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundFragment.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundFragment.cc	(revision 1007)
@@ -24,9 +24,7 @@
 // ********************************************************************
 //
-// $Id: G4VPreCompoundFragment.cc,v 1.12 2009/02/10 16:01:37 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
 //
-// J. M. Quesada (August 2008).  
-// Based  on previous work by V. Lara
+//J. M. Quesada (August 2008).  
+//Based  on previous work by V. Lara
 //
  
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundIon.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundIon.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundIon.cc	(revision 1007)
@@ -26,5 +26,5 @@
 //
 // $Id: G4VPreCompoundIon.cc,v 1.9 2008/09/22 10:18:36 ahoward Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+// GEANT4 tag $Name: geant4-09-02 $
 //
 // by V. Lara
Index: trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundNucleon.cc
===================================================================
--- trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundNucleon.cc	(revision 968)
+++ trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4VPreCompoundNucleon.cc	(revision 1007)
@@ -27,5 +27,5 @@
 
 // $Id: G4VPreCompoundNucleon.cc,v 1.9 2008/09/22 10:18:36 ahoward Exp $
-// GEANT4 tag $Name: geant4-09-02-ref-02 $
+// GEANT4 tag $Name: geant4-09-02 $
 
 //
