#include <G4ScreeningMottCrossSection.hh>

Collaboration diagram for G4ScreeningMottCrossSection:

Public Member Functions
	G4ScreeningMottCrossSection ()

virtual	~G4ScreeningMottCrossSection ()

void	Initialise (const G4ParticleDefinition *, G4double cosThetaLim)

G4double	GetScreeningAngle ()

void	SetScreeningCoefficient ()

void	SetupParticle (const G4ParticleDefinition *)

void	SetupKinematic (G4double kinEnergy, G4double Z)

G4double	NuclearCrossSection ()

G4ThreeVector	GetNewDirection ()

G4double	GetMom2CM () const

G4double	GetMom2Lab () const

G4double	GetTrec () const

G4double	GetScreeningCoefficient () const

G4double	GetTotalCross () const

G4double	McFcorrection (G4double)

G4double	RatioMottRutherford (G4double)

G4double	FormFactor2ExpHof (G4double)

G4double	GetScatteringAngle ()

G4double	AngleDistribution (G4double)

Private Member Functions
G4ScreeningMottCrossSection &	operator= (const G4ScreeningMottCrossSection &right)

	G4ScreeningMottCrossSection (const G4ScreeningMottCrossSection &)

Private Attributes
G4NistManager *	fNistManager

G4Pow *	fG4pow

G4MottCoefficients *	mottcoeff

G4double	TotalCross

G4double	cosThetaMin

G4double	cosThetaMax

G4double	cosTetMinNuc

G4double	cosTetMaxNuc

G4double	ecut

G4double	etag

const G4ParticleDefinition *	particle

G4double	spin

G4double	mass

G4double	tkinLab

G4double	momLab2

G4double	invbetaLab2

G4double	mu_rel

G4double	tkin

G4double	mom2

G4double	invbeta2

G4double	beta

G4double	gamma

G4double	targetZ

G4double	targetA

G4double	targetMass

G4double	Trec

G4double	As

G4double	coeffb [5][6]

G4double	alpha

G4double	htc2

G4double	e2

G4double	cross [DIM]

Detailed Description

Definition at line 88 of file G4ScreeningMottCrossSection.hh.

Constructor & Destructor Documentation

◆ G4ScreeningMottCrossSection() [1/2]

G4ScreeningMottCrossSection::G4ScreeningMottCrossSection ( )

Definition at line 86 of file G4ScreeningMottCrossSection.cc.

                                                         :
    cosThetaMin(1.0),
    cosThetaMax(-1.0),
    alpha(fine_structure_const),
    htc2(hbarc_squared),
    e2(electron_mass_c2*classic_electr_radius) 
 {
   TotalCross=0;
 
   fNistManager = G4NistManager::Instance();
   fG4pow = G4Pow::GetInstance();
   particle=0;
 
   spin = mass = mu_rel=0;
   tkinLab = momLab2 = invbetaLab2=0;
   tkin = mom2 = invbeta2=beta=gamma=0;
 
   Trec=targetZ = targetMass = As =0;
   etag = ecut = 0.0;
 
   targetA = 0;
 
   cosTetMinNuc=0;
   cosTetMaxNuc=0;
 
   for(G4int i=0 ; i<5;  ++i){
     for(G4int j=0; j< 6; ++j){
       coeffb[i][j]=0;
     }
   }
   if(dangle[0] == 0.0) {
     for(G4int i=0; i<DIM; ++i){
       cross[i]=0;
     }
   
  
   }
 
   mottcoeff = new G4MottCoefficients();
 }

Here is the call graph for this function:

◆ ~G4ScreeningMottCrossSection()

G4ScreeningMottCrossSection::~G4ScreeningMottCrossSection ( )

virtual

Definition at line 129 of file G4ScreeningMottCrossSection.cc.

 {
   delete mottcoeff;
 }

◆ G4ScreeningMottCrossSection() [2/2]

G4ScreeningMottCrossSection::G4ScreeningMottCrossSection ( const G4ScreeningMottCrossSection & )

private

Member Function Documentation

◆ AngleDistribution()

G4double G4ScreeningMottCrossSection::AngleDistribution ( G4double anglein )

Definition at line 329 of file G4ScreeningMottCrossSection.cc.

 {
   G4double total=TotalCross ;
   G4double fatt= e2*targetZ/(mu_rel*gamma*beta*beta);
   G4double fatt2=fatt*fatt;
   total/=fatt2;
 
   G4double R=0;
   if (coeffb[0][0]!=0){
     //   cout<<" Mott....targetZ "<< targetZ<<endl;      
     R=RatioMottRutherford(anglein);
   } else if (coeffb[0][0]==0){
     // cout<<" McF.... targetZ "<< targetZ<<endl;
     R=McFcorrection(anglein);
   }
 
   G4double y=twopi*sin(anglein)*R*FormFactor2ExpHof(anglein)/
     ( (2*As+2.*fG4pow->powN(sin(anglein*0.5),2))*(2.*As+2.*fG4pow->powN(sin(anglein*0.5),2)) );
   //G4cout<<"AngleDistribution Total Cross: "<<TotalCross<<" y: "<<y<<G4endl;//TO REMOVE
   return y/total;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ FormFactor2ExpHof()

G4double G4ScreeningMottCrossSection::FormFactor2ExpHof ( G4double angles )

Definition at line 233 of file G4ScreeningMottCrossSection.cc.

 {
   G4double M=targetMass; 
   G4double E=tkinLab;
   G4double Etot=E+mass;
   G4double Tmax=2.*M*E*(E+2.*mass)/(mass*mass+M*M+2.*M*Etot);
   G4double T=Tmax*fG4pow->powN(sin(angles*0.5), 2);
   G4double q2=T*(T+2.*M);
   q2/=htc2;//1/cm2
   G4double RN=1.27e-13*G4Exp(G4Log(targetA)*0.27)*cm;
   G4double xN= (RN*RN*q2);
   G4double den=(1.+xN/12.);
   G4double FN=1./(den*den);
   G4double form2=(FN*FN);
 
   return form2;
 
   //cout<<"..................... form2 "<< form2<<endl;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ GetMom2CM()

G4double G4ScreeningMottCrossSection::GetMom2CM ( ) const

inline

Definition at line 200 of file G4ScreeningMottCrossSection.hh.

 {
     return mom2;
 }

◆ GetMom2Lab()

G4double G4ScreeningMottCrossSection::GetMom2Lab ( ) const

inline

Definition at line 207 of file G4ScreeningMottCrossSection.hh.

 {
         return momLab2;
 }

◆ GetNewDirection()

G4ThreeVector G4ScreeningMottCrossSection::GetNewDirection ( )

Definition at line 381 of file G4ScreeningMottCrossSection.cc.

                                                           {
 
   G4ThreeVector dir(0.0,0.0,1.0);
     
   G4double z1=GetScatteringAngle();
   
   G4double sint = sin(z1);
   G4double cost = sqrt(1.0 - sint*sint);
   G4double phi  = twopi* G4UniformRand();
   G4double dirx = sint*cos(phi);
   G4double diry = sint*sin(phi);
   G4double dirz = cost;
 
   //.......set Trc
   G4double etot=tkinLab+mass;
   G4double mass2=targetMass;
   Trec=(1.0 - cost)* mass2*(etot*etot - mass*mass )/
     (mass*mass + mass2*mass2+ 2.*mass2*etot);
        
   dir.set(dirx,diry,dirz);
 
   return dir;
 }

Here is the call graph for this function:

◆ GetScatteringAngle()

G4double G4ScreeningMottCrossSection::GetScatteringAngle ( )

Definition at line 353 of file G4ScreeningMottCrossSection.cc.

 {
   //    cout<<"................ tkinLab  "<< G4BestUnit(tkinLab,"Energy") << " anglemin=  "<<anglemin<<endl;
   //cout<<"anglemax=  "<<anglemax<<endl;
   G4double r =G4UniformRand();
   
   G4double scattangle=0;
   G4double y=0;
   G4double dy=0;
   G4double area=0;
 
   for(G4int i=0; i<DIM; ++i){
     
     y+=AngleDistribution(tet[i])*dangle[i];
     dy= y-area ;
     area=y;
 
     if(r >=y-dy && r<y ){   
       scattangle= angle[i] +G4UniformRand()*dangle[i];
       //G4cout<<"Energy: "<<tkinLab/MeV<<" Random: "<<r<<" y: "<<y<<" Scattangle: "<<scattangle<<" Angle: "<<angle[i]<<" dAngle: "<<dangle[i]<<G4endl; //TO REMOVE
       break;
     }           
   }
   return scattangle;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ GetScreeningAngle()

G4double G4ScreeningMottCrossSection::GetScreeningAngle ( )

Definition at line 162 of file G4ScreeningMottCrossSection.cc.

 {
   SetScreeningCoefficient();
 
   G4double screenangle=2.*asin(sqrt(As));
   //    cout<<"  screenangle  "<<  screenangle <<endl;
   if(screenangle>=pi) screenangle=pi;
     
   return screenangle;
 }

Here is the call graph for this function:

◆ GetScreeningCoefficient()

G4double G4ScreeningMottCrossSection::GetScreeningCoefficient ( ) const

inline

Definition at line 222 of file G4ScreeningMottCrossSection.hh.

 {
     return As;
 }

◆ GetTotalCross()

G4double G4ScreeningMottCrossSection::GetTotalCross ( ) const

inline

Definition at line 230 of file G4ScreeningMottCrossSection.hh.

 {
     return TotalCross;
 }

◆ GetTrec()

G4double G4ScreeningMottCrossSection::GetTrec ( ) const

inline

Definition at line 215 of file G4ScreeningMottCrossSection.hh.

 {
         return Trec;
 }

◆ Initialise()

void G4ScreeningMottCrossSection::Initialise	(	const G4ParticleDefinition *	p,
		G4double	cosThetaLim
	)

Definition at line 135 of file G4ScreeningMottCrossSection.cc.

 {
   SetupParticle(p);
   tkin = targetZ = mom2 = DBL_MIN;
   ecut = etag = DBL_MAX;
   particle = p;
   cosThetaMin = CosThetaLim; 
 
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ McFcorrection()

G4double G4ScreeningMottCrossSection::McFcorrection ( G4double angles )

Definition at line 255 of file G4ScreeningMottCrossSection.cc.

 {
   G4double  beta2=1./invbeta2;
   G4double sintmezzi=std::sin(angles/2.);
   G4double sin2tmezzi = sintmezzi*sintmezzi;
   G4double R=1.-beta2*sin2tmezzi + targetZ*alpha*beta*pi*sintmezzi*(1.-sintmezzi);
   return R;
 }

Here is the caller graph for this function:

◆ NuclearCrossSection()

G4double G4ScreeningMottCrossSection::NuclearCrossSection ( )

Definition at line 290 of file G4ScreeningMottCrossSection.cc.

 {
   if(cosTetMaxNuc >= cosTetMinNuc) return 0.0;
 
   TotalCross=0;
 
   for(G4int i=0; i<DIM; ++i){
     G4double R=0;
     G4double F2=FormFactor2ExpHof(tet[i]);
             
     if (coeffb[0][0]!=0){
       //cout<<" Mott....targetZ "<< targetZ<<endl;  
       R=RatioMottRutherford(tet[i]);
     } else if (coeffb[0][0]==0){
       // cout<<" McF.... targetZ "<< targetZ<<endl;
       R=McFcorrection(tet[i]);
     }
 
     //cout<<"----------------- R "<<R<<" F2 "<<F2<<endl;
     //                cout<<"angle "<<tet[i] << " F2 "<<F2<<endl;
 
     G4double e4=e2*e2;
     G4double den=2.*As+2.*fG4pow->powN(sin(tet[i]*0.5),2);
     G4double func=1./(den*den);
 
     G4double fatt= targetZ/(mu_rel*gamma*beta*beta);
     G4double sigma=e4*fatt*fatt*func;
     cross[i]=F2*R*sigma;
     G4double pi2sintet=twopi*sin(tet[i]);
 
     TotalCross+=pi2sintet*cross[i]*dangle[i];
   }//end integral
 
   //cout<< "ok ......... TotalCross "<<TotalCross<<endl;
   return TotalCross;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ operator=()

G4ScreeningMottCrossSection& G4ScreeningMottCrossSection::operator= ( const G4ScreeningMottCrossSection & right )

private

◆ RatioMottRutherford()

G4double G4ScreeningMottCrossSection::RatioMottRutherford ( G4double angles )

Definition at line 264 of file G4ScreeningMottCrossSection.cc.

 {
   G4double R=0;
   G4double fcost=std::sqrt((1. -cos(angles)));
   G4double a[5];
   const G4double shift=0.7181228;
   G4double beta0= beta -shift;
 
   for(G4int j=0 ;j<=4;j++){
     a[j]=0;
   }
 
   for(G4int j=0 ;j<=4;j++){
     for(G4int k=0;k<=5;k++ ){  
       a[j]+=coeffb[j][k]*fG4pow->powN(beta0,k);
     }
   }
 
   for(G4int j=0; j<=4; ++j){
     R+=a[j]*fG4pow->powN(fcost,j);
   }
   return R;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ SetScreeningCoefficient()

void G4ScreeningMottCrossSection::SetScreeningCoefficient ( )

Definition at line 146 of file G4ScreeningMottCrossSection.cc.

 {
   G4double alpha2=alpha*alpha;
   //Bohr radius 
   G4double a0=  Bohr_radius  ;//0.529e-8*cm;
   //Thomas-Fermi screening length
   G4double aU=0.88534*a0/fG4pow->Z13(targetZ);
   G4double twoR2=aU*aU;
 
   G4double factor= 1.13 + 3.76*targetZ*targetZ*invbeta2*alpha2;
   As=0.25*(htc2)/(twoR2*mom2)*factor;
   //cout<<"0k .........................As  "<<As<<endl;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ SetupKinematic()

void G4ScreeningMottCrossSection::SetupKinematic	(	G4double	kinEnergy,
		G4double	Z
	)

Definition at line 174 of file G4ScreeningMottCrossSection.cc.

 {
   //...Target
   G4int iz = G4int(Z);
   G4double A = fNistManager->GetAtomicMassAmu(iz);
   G4int ia = G4int(A);
   G4double mass2 = G4NucleiProperties::GetNuclearMass(ia, iz);
 
   targetZ = Z;
   targetA = fNistManager->GetAtomicMassAmu(iz);
   targetMass= mass2;
 
   mottcoeff->SetMottCoeff(targetZ, coeffb);
 
   //cout<<"......... targetA "<< targetA <<endl;
   //cout<<"......... targetMass "<< targetMass/MeV <<endl;
 
   // incident particle lab
   tkinLab = ekin;
   momLab2 = tkinLab*(tkinLab + 2.0*mass);
   invbetaLab2 = 1.0 +  mass*mass/momLab2;
 
   G4double etot = tkinLab + mass;
   G4double ptot = sqrt(momLab2);
   G4double m12  = mass*mass;
                 
   // relativistic reduced mass from publucation
   // A.P. Martynenko, R.N. Faustov, Teoret. mat. Fiz. 64 (1985) 179
         
   //incident particle & target nucleus
   G4double Ecm=sqrt(m12 + mass2*mass2 + 2.0*etot*mass2);
   mu_rel=mass*mass2/Ecm;
   G4double momCM= ptot*mass2/Ecm;
   // relative system
   mom2 = momCM*momCM;
   invbeta2 = 1.0 +  mu_rel*mu_rel/mom2;
   tkin = momCM*sqrt(invbeta2) - mu_rel;//Ekin of mu_rel
   G4double  beta2=1./invbeta2;
   beta=std::sqrt(beta2) ;
   G4double gamma2= 1./(1.-beta2);
   gamma=std::sqrt(gamma2);
 
   //.........................................................
 
   
   SetScreeningCoefficient();
 
   //Integration Angles definition
 
   cosTetMinNuc =cosThetaMin;
   cosTetMaxNuc =cosThetaMax;
 
   for(G4int i=0; i<DIM;  ++i ){
     cross[i]=0;
   }
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ SetupParticle()

void G4ScreeningMottCrossSection::SetupParticle ( const G4ParticleDefinition * p )

inline

Definition at line 189 of file G4ScreeningMottCrossSection.hh.

 {
     particle = p;
     mass = particle->GetPDGMass();
     spin = particle->GetPDGSpin();
     if(0.0 != spin) { spin = 0.5; }
     tkin = 0.0;
 }