G4ScreeningMottCrossSection Class Reference

#include <G4ScreeningMottCrossSection.hh>

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)

Detailed Description

Definition at line 86 of file G4ScreeningMottCrossSection.hh.

Constructor & Destructor Documentation

G4ScreeningMottCrossSection::G4ScreeningMottCrossSection

(

)

Definition at line 83 of file G4ScreeningMottCrossSection.cc.

References G4NistManager::Instance().

                                                        :
   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();
        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;
        }
    }

    mottcoeff = new G4MottCoefficients();
}

G4ScreeningMottCrossSection::~G4ScreeningMottCrossSection ( )

virtual

Definition at line 119 of file G4ScreeningMottCrossSection.cc.

{
  delete mottcoeff;
}

Member Function Documentation

G4double G4ScreeningMottCrossSection::AngleDistribution

(

G4double

anglein

)

Definition at line 394 of file G4ScreeningMottCrossSection.cc.

References FormFactor2ExpHof(), McFcorrection(), python.hepunit::pi, RatioMottRutherford(), and G4INCL::CrossSections::total().

Referenced by GetScatteringAngle().

                                                                       {


    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=2.*pi*sin(anglein)*R*FormFactor2ExpHof(anglein)/
                ((2*As+2.*pow(sin(anglein/2.),2.))*(2.*As+2.*pow(sin(anglein/2.),2.) ));

return y/total;

}

G4double G4ScreeningMottCrossSection::FormFactor2ExpHof

(

G4double

angle

)

Definition at line 238 of file G4ScreeningMottCrossSection.cc.

References python.hepunit::cm, and form2().

Referenced by AngleDistribution(), and NuclearCrossSection().

                                                                     {


       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*pow(sin(angle/2.),2.);
       G4double q2=T*(T+2.*M);
             q2/=htc2;//1/cm2
       G4double RN=1.27e-13*pow(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;
}

G4double G4ScreeningMottCrossSection::GetMom2CM ( ) const

inline

Definition at line 195 of file G4ScreeningMottCrossSection.hh.

{
    return mom2;
}

G4double G4ScreeningMottCrossSection::GetMom2Lab ( ) const

inline

Definition at line 202 of file G4ScreeningMottCrossSection.hh.

Referenced by G4eSingleCoulombScatteringModel::SampleSecondaries().

{
        return momLab2;
}

G4ThreeVector G4ScreeningMottCrossSection::GetNewDirection

(

)

Definition at line 497 of file G4ScreeningMottCrossSection.cc.

References G4UniformRand, GetScatteringAngle(), CLHEP::Hep3Vector::set(), and python.hepunit::twopi.

Referenced by G4eSingleCoulombScatteringModel::SampleSecondaries().

                                                          {

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;

}

G4double G4ScreeningMottCrossSection::GetScatteringAngle

(

)

Definition at line 423 of file G4ScreeningMottCrossSection.cc.

References AngleDistribution(), G4UniformRand, and GetScreeningAngle().

Referenced by GetNewDirection().

{
  

//cout<<" G4ScreeningMottCrossSection::SampleCosineTheta ............."<<endl;  

    if(cosTetMaxNuc >= cosTetMinNuc) return 0.0;

    G4double anglemin=std::acos(cosTetMinNuc);  
    G4double anglemax= std::acos(cosTetMaxNuc);

        static const G4double limit = 1.e-9;
        if(anglemin < limit) {
      anglemin = GetScreeningAngle()/10.;
          if(anglemin < limit) { anglemin = limit; }
        }

//  cout<<"................ tkinLab  "<< G4BestUnit(tkinLab,"Energy") << " anglemin=  "<<anglemin<<endl;
    //cout<<"anglemax=  "<<anglemax<<endl;
    G4double r =G4UniformRand();

    G4double loganglemin=log10(anglemin);
        G4double loganglemax=log10(anglemax);
        G4double logdangle=0.01;

        G4int bins=(G4int)((loganglemax-loganglemin)/logdangle);

        std::vector<G4double> angle;
        std::vector<G4double> tet;
        std::vector<G4double> dangle;

              for(G4int i=0; i<=bins; i++ ){
                     tet.push_back(0);
                     dangle.push_back(0);
                     angle.push_back(pow(10.,loganglemin+logdangle*i));
                        }


                G4int  dim = tet.size();
        G4double scattangle=0;
        G4double y=0;
        G4double dy=0;
        G4double area=0;

                for(G4int i=0; i<dim;i++){

                        if(i!=dim-1){
                        dangle[i]=(angle[i+1]-angle[i]);
                        tet[i]=(angle[i+1]+angle[i])/2.;
                        }else if(i==dim-1){
                                break;
            }

            y+=AngleDistribution(tet[i])*dangle[i];
            dy= y-area ;
            area=y;

            if(r >=y-dy && r<=y+dy ){   
            scattangle= angle[i] +G4UniformRand()*dangle[i];
            //cout<<"y "<<y <<" r  "<< r  << " y/r "<<y/r << endl;
            break;

            }   
        
        }


  return scattangle;


}

G4double G4ScreeningMottCrossSection::GetScreeningAngle

(

)

Definition at line 158 of file G4ScreeningMottCrossSection.cc.

References python.hepunit::pi, and SetScreeningCoefficient().

Referenced by GetScatteringAngle(), NuclearCrossSection(), and SetupKinematic().

                                                       {

                SetScreeningCoefficient();

        //cout<<" As "<<As<<endl;
                if(As < 0.0) { As = 0.0; }
                else if(As > 1.0) { As = 1.0; }

                G4double screenangle=2.*asin(sqrt(As));

    //  cout<<"  screenangle  "<<  screenangle <<endl;
    
        if(screenangle>=pi) screenangle=pi;

    
return screenangle;

}

G4double G4ScreeningMottCrossSection::GetScreeningCoefficient ( ) const

inline

Definition at line 217 of file G4ScreeningMottCrossSection.hh.

{
    return As;
}

G4double G4ScreeningMottCrossSection::GetTotalCross ( ) const

inline

Definition at line 225 of file G4ScreeningMottCrossSection.hh.

Referenced by G4eSingleCoulombScatteringModel::SampleSecondaries().

{
    return TotalCross;
}

G4double G4ScreeningMottCrossSection::GetTrec ( ) const

inline

Definition at line 210 of file G4ScreeningMottCrossSection.hh.

Referenced by G4eSingleCoulombScatteringModel::SampleSecondaries().

{
        return Trec;
}

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

Definition at line 125 of file G4ScreeningMottCrossSection.cc.

References DBL_MAX, DBL_MIN, and SetupParticle().

Referenced by G4eSingleCoulombScatteringModel::Initialise().

{

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

}

G4double G4ScreeningMottCrossSection::McFcorrection

(

G4double

angle

)

Definition at line 262 of file G4ScreeningMottCrossSection.cc.

References python.hepunit::pi.

Referenced by AngleDistribution(), and NuclearCrossSection().

                                                                  {


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


return R;
}

G4double G4ScreeningMottCrossSection::NuclearCrossSection

(

)

Definition at line 304 of file G4ScreeningMottCrossSection.cc.

References FormFactor2ExpHof(), GetScreeningAngle(), McFcorrection(), python.hepunit::pi, and RatioMottRutherford().

Referenced by G4eSingleCoulombScatteringModel::ComputeCrossSectionPerAtom().

{
     if(cosTetMaxNuc >= cosTetMinNuc) return 0.0;

        TotalCross=0;

        G4double anglemin =std::acos(cosTetMinNuc);
        G4double anglemax =std::acos(cosTetMaxNuc); 

        static const G4double limit = 1.e-9;
        if(anglemin < limit) {
      anglemin = GetScreeningAngle()/10.;
          if(anglemin < limit) { anglemin = limit; }
        }

    //cout<<" anglemin  "<< anglemin <<endl;

        G4double loganglemin=log10(anglemin);
        G4double loganglemax=log10(anglemax);
        G4double logdangle=0.01;

        G4int bins=(G4int)((loganglemax-loganglemin)/logdangle);


        vector<G4double> angle;
        vector<G4double> tet;
        vector<G4double> dangle;
        vector<G4double> cross;

              for(G4int i=0; i<=bins; i++ ){
                     tet.push_back(0);
                     dangle.push_back(0);
                     angle.push_back(pow(10.,loganglemin+logdangle*i));
                     cross.push_back(0);
                }


        G4int  dim = tet.size();
        //cout<<"dim--- "<<dim<<endl;


                for(G4int i=0; i<dim;i++){

                        if(i!=dim-1){
                        dangle[i]=(angle[i+1]-angle[i]);
                        tet[i]=(angle[i+1]+angle[i])/2.;
                    }else if(i==dim-1){
                                break;
            }


                    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.*pow(sin(tet[i]/2.),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=2.*pi*sin(tet[i]);

                    TotalCross+=pi2sintet*cross[i]*dangle[i];
                }//end integral


//cout<< "ok ......... TotalCross "<<TotalCross<<endl;


return TotalCross;

}

G4double G4ScreeningMottCrossSection::RatioMottRutherford

(

G4double

angle

)

Definition at line 274 of file G4ScreeningMottCrossSection.cc.

References test::a.

Referenced by AngleDistribution(), and NuclearCrossSection().

                                                                       {


        G4double R=0;
        G4double fcost=std::sqrt((1. -cos(angle)));
        G4double a[5];
        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]*pow(beta0,k);
                    }
            }

    for(G4int j=0 ;j<=4 ;j++){
                R+=a[j]* pow(fcost,j);
                }



return R;

}

void G4ScreeningMottCrossSection::SetScreeningCoefficient

(

)

Definition at line 137 of file G4ScreeningMottCrossSection.cc.

References python.hepunit::Bohr_radius.

Referenced by GetScreeningAngle().

                                                         {

        G4double alpha2=alpha*alpha;
        //Bohr radius 
        G4double a0=  Bohr_radius  ;//0.529e-8*cm;
        //Thomas-Fermi screening length
        G4double aU=0.88534*a0/pow(targetZ,1./3.);
        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;

}

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

Definition at line 178 of file G4ScreeningMottCrossSection.cc.

References G4NistManager::GetAtomicMassAmu(), G4NucleiProperties::GetNuclearMass(), GetScreeningAngle(), iz, G4INCL::Math::min(), and G4MottCoefficients::SetMottCoeff().

Referenced by G4eSingleCoulombScatteringModel::ComputeCrossSectionPerAtom().

{

        //...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);



               //.........................................................


        G4double screenangle=GetScreeningAngle()/10.;
        //cout<<" screenangle [rad] "<<screenangle/rad <<endl;

                cosTetMinNuc =min( cosThetaMin ,cos(screenangle));
                cosTetMaxNuc =cosThetaMax;
    
//cout<<"ok..............mu_rel "<<mu_rel<<endl;

}

void G4ScreeningMottCrossSection::SetupParticle ( const G4ParticleDefinition *  p )

inline

Definition at line 184 of file G4ScreeningMottCrossSection.hh.

References G4ParticleDefinition::GetPDGMass(), and G4ParticleDefinition::GetPDGSpin().

Referenced by Initialise().

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

---

The documentation for this class was generated from the following files:

• G4ScreeningMottCrossSection.hh
• G4ScreeningMottCrossSection.cc