Changeset 1340 for trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc

Timestamp:

Nov 5, 2010, 3:45:55 PM (14 years ago)

Author:

garnier

Message:

update ti head

File:

• trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc (modified) (2 diffs)

trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundNucleon.cc

@@ -24 +24 @@
 // ********************************************************************
 //
-//
-// $Id: G4PreCompoundNucleon.cc,v 1.13 2009/02/11 18:06:00 vnivanch Exp $
-// GEANT4 tag $Name: geant4-09-04-beta-01 $
+// $Id: G4PreCompoundNucleon.cc,v 1.14 2010/08/28 15:16:55 vnivanch Exp $
+// GEANT4 tag $Name: geant4-09-03-ref-09 $
 //
 // -------------------------------------------------------------------
@@ -39 +38 @@
 // Modified:  
 // 10.02.2009 J. M. Quesada cleanup
+// 20.08.2010 V.Ivanchenko added G4Pow and G4PreCompoundParameters pointers
+//                         use int Z and A and cleanup
 //
 
 #include "G4PreCompoundNucleon.hh"
-#include "G4PreCompoundParameters.hh"
 
-G4bool G4PreCompoundNucleon::IsItPossible(const G4Fragment& aFragment) 
+G4PreCompoundNucleon::
+G4PreCompoundNucleon(const G4ParticleDefinition* part,
+                      G4VCoulombBarrier* aCoulombBarrier)
+  : G4PreCompoundFragment(part,aCoulombBarrier)
 {
-  G4int pplus = aFragment.GetNumberOfCharged();   
-  G4int pneut = aFragment.GetNumberOfParticles()-pplus;
-  return (pneut >= (GetA()-GetZ()) && pplus >= GetZ());
+  fact = 2*CLHEP::millibarn/(CLHEP::pi2*CLHEP::hbarc*CLHEP::hbarc*CLHEP::hbarc);
 }
 
+G4PreCompoundNucleon::~G4PreCompoundNucleon()
+{}
+
 G4double G4PreCompoundNucleon::
-ProbabilityDistributionFunction(const G4double eKin, 
+ProbabilityDistributionFunction(G4double eKin, 
                                 const G4Fragment& aFragment)
 {
-  if ( !IsItPossible(aFragment) ) return 0.0;
+  if ( !IsItPossible(aFragment) ) { return 0.0; }
 
   G4double U = aFragment.GetExcitationEnergy();
-  G4double P = aFragment.GetNumberOfParticles();
-  G4double H = aFragment.GetNumberOfHoles();
-  G4double N = P + H;
+  G4int P = aFragment.GetNumberOfParticles();
+  G4int H = aFragment.GetNumberOfHoles();
+  G4int N = P + H;
+
+  G4double g0 = (6.0/pi2)*aFragment.GetA_asInt()*theParameters->GetLevelDensity();
+  G4double g1 = (6.0/pi2)*GetRestA()*theParameters->GetLevelDensity();
   
-  G4double g0 = (6.0/pi2)*aFragment.GetA() * 
-    G4PreCompoundParameters::GetAddress()->GetLevelDensity();
- 
-  G4double g1 = (6.0/pi2)*GetRestA() * 
-    G4PreCompoundParameters::GetAddress()->GetLevelDensity();
+  G4double A0 = G4double(P*P+H*H+P-3*H)/(4.0*g0);
+  G4double A1 = (A0 - 0.5*P)/g1;  
 
-  G4double A0 = ((P*P+H*H+P-H)/4.0 - H/2.0)/g0;
+  G4double E0 = U - A0;
+  if (E0 <= 0.0) { return 0.0; }
 
-  G4double A1 = (A0 - P/2.0)/g1;
+  G4double E1 = U - eKin - GetBindingEnergy() - A1;
+  if (E1 <= 0.0) { return 0.0; }
+
+  G4double rj = GetRj(P, aFragment.GetNumberOfCharged());
+  G4double xs = CrossSection(eKin);
+
+  if (rj <0.0 || xs < 0.0) {  
+    std::ostringstream errOs;
+    G4cout<<"WARNING:  NEGATIVE VALUES "<<G4endl;     
+    errOs << "Rj=" << rj <<"  xsec("
+          <<eKin/MeV<<" MeV)= "<< xs <<"  A= "<<GetA()<<"  Z= "<<GetZ()
+          <<G4endl;
+    throw G4HadronicException(__FILE__, __LINE__, errOs.str());
+  }
+
+  G4double Probability = fact * GetReducedMass() * rj * xs * eKin * P * (N-1) 
+    * g4pow->powN(g1*E1/(g0*E0),N-2) * g1 / (E0*g0*g0);
   
-  G4double E0 = std::max(0.0,U - A0);
-  if (E0 == 0.0) return 0.0;
-  G4double E1 = std::max(0.0,U - eKin - GetBindingEnergy() - A1);
-  if (E1 == 0.0) return 0.0;
-
-
+  /*
   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);
-
-  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());
-  }
+    * eKin*CrossSection(eKin)*millibarn * P*(N-1.0) * 
+   std::pow(g1*E1/(g0*E0),N-2.0)/E0 * g1/(g0*g0);
+  */
 
   return Probability;