PrimaryGeneratorAction2 Class Reference

#include <PrimaryGeneratorAction2.hh>

Collaboration diagram for PrimaryGeneratorAction2:

Public Member Functions
	PrimaryGeneratorAction2 (G4ParticleGun *)
	~PrimaryGeneratorAction2 ()
void	GeneratePrimaries (G4Event *)
G4double	RejectAccept ()
G4double	InverseCumul ()

Private Member Functions
void	InitFunction ()

Private Attributes
G4ParticleGun *	particleGun
G4int	nPoints
std::vector< G4double >	x
std::vector< G4double >	f
std::vector< G4double >	a
std::vector< G4double >	Fc
G4double	fMax

Detailed Description

Definition at line 47 of file PrimaryGeneratorAction2.hh.

Constructor & Destructor Documentation

PrimaryGeneratorAction2()

PrimaryGeneratorAction2::PrimaryGeneratorAction2

(

G4ParticleGun *

gun

)

Definition at line 48 of file PrimaryGeneratorAction2.cc.

: particleGun(gun)
{    
  // energy distribution
  //
  InitFunction();
}

Here is the call graph for this function:

~PrimaryGeneratorAction2()

PrimaryGeneratorAction2::~PrimaryGeneratorAction2

(

)

Definition at line 58 of file PrimaryGeneratorAction2.cc.

{ }

Member Function Documentation

GeneratePrimaries()

void PrimaryGeneratorAction2::GeneratePrimaries

(

G4Event *

anEvent

)

Definition at line 63 of file PrimaryGeneratorAction2.cc.

{
  //cosAlpha uniform in [cos(0), cos(pi)]
  G4double cosAlpha = 1. - 2*G4UniformRand();
  G4double sinAlpha = std::sqrt(1. - cosAlpha*cosAlpha);
  G4double psi      = twopi*G4UniformRand();  //psi uniform in [0, 2*pi]  
  G4ThreeVector dir(sinAlpha*std::cos(psi),sinAlpha*std::sin(psi),cosAlpha);

  particleGun->SetParticleMomentumDirection(dir);
  
  //set energy from a tabulated distribution
  //
  //G4double energy = RejectAccept();
  G4double energy = InverseCumul();  
  particleGun->SetParticleEnergy(energy);    

  //create vertex
  //   
  particleGun->GeneratePrimaryVertex(anEvent);
}

Here is the call graph for this function:

InitFunction()

void PrimaryGeneratorAction2::InitFunction ( )

private

Definition at line 85 of file PrimaryGeneratorAction2.cc.

{
  // tabulated function 
  // f is assumed positive, linear per segment, continuous
  //
  nPoints = 16;
  const G4double xx[] = 
    { 37*keV, 39*keV, 45*keV,  51*keV,  57*keV,  69*keV,  71*keV,  75*keV, 
      83*keV, 91*keV, 97*keV, 107*keV, 125*keV, 145*keV, 159*keV, 160*keV }; 
      
  const G4double ff[] =
    { 0.000,  0.077,  0.380,  2.044, 5.535, 15.077, 12.443, 14.766,
     17.644, 18.518, 17.772, 14.776, 8.372,  3.217,  0.194,  0.000 };
  
  //copy arrays in std::vector and compute fMax
  //
  x.resize(nPoints); f.resize(nPoints);
  fMax = 0.;
  for (G4int j=0; j<nPoints; j++) {
    x[j] = xx[j]; f[j] = ff[j];
    if (fMax < f[j]) fMax = f[j];
  };
     
  //compute slopes
  //
  a.resize(nPoints);
  for (G4int j=0; j<nPoints-1; j++) { 
    a[j] = (f[j+1] - f[j])/(x[j+1] - x[j]);
  };
  
  //compute cumulative function
  //
  Fc.resize(nPoints);  
  Fc[0] = 0.;
  for (G4int j=1; j<nPoints; j++) {
    Fc[j] = Fc[j-1] + 0.5*(f[j] + f[j-1])*(x[j] - x[j-1]);
  };     
}

Here is the caller graph for this function:

InverseCumul()

G4double PrimaryGeneratorAction2::InverseCumul

(

)

Definition at line 148 of file PrimaryGeneratorAction2.cc.

{
  // tabulated function
  // f is assumed positive, linear per segment, continuous 
  // --> cumulative function is second order polynomial
  
  //choose y randomly
  G4double y_rndm = G4UniformRand()*Fc[nPoints-1];
  //find bin
  G4int j = nPoints-2;
  while ((Fc[j] > y_rndm) && (j > 0)) j--;
  //y_rndm --> x_rndm :  Fc(x) is second order polynomial
  G4double x_rndm = x[j];
  G4double aa = a[j];
  if (aa != 0.) {
    G4double b = f[j]/aa, c = 2*(y_rndm - Fc[j])/aa;
    G4double delta = b*b + c;
    G4int sign = 1; if (aa < 0.) sign = -1;
    x_rndm += sign*std::sqrt(delta) - b;    
  } else if (f[j] > 0.) {
    x_rndm += (y_rndm - Fc[j])/f[j];
  };
  return x_rndm;
}

Here is the call graph for this function:

Here is the caller graph for this function:

RejectAccept()

G4double PrimaryGeneratorAction2::RejectAccept

(

)

Definition at line 126 of file PrimaryGeneratorAction2.cc.

{
  // tabulated function 
  // f is assumed positive, linear per segment, continuous
  //  
  G4double x_rndm = 0., y_rndm = 0., f_inter = -1.;
  
  while (y_rndm > f_inter) {
    //choose a point randomly
    x_rndm = x[0] + G4UniformRand()*(x[nPoints-1] - x[0]);
    y_rndm = G4UniformRand()*fMax;
    //find bin
    G4int j = nPoints-2;
    while ((x[j] > x_rndm) && (j > 0)) j--;
    //compute f(x_rndm) by linear interpolation
    f_inter = f[j] + a[j]*(x_rndm - x[j]);
  };
  return x_rndm;
}

Member Data Documentation

a

std::vector<G4double> PrimaryGeneratorAction2::a

private

Definition at line 66 of file PrimaryGeneratorAction2.hh.

f

std::vector<G4double> PrimaryGeneratorAction2::f

private

Definition at line 65 of file PrimaryGeneratorAction2.hh.

Fc

std::vector<G4double> PrimaryGeneratorAction2::Fc

private

Definition at line 67 of file PrimaryGeneratorAction2.hh.

fMax

G4double PrimaryGeneratorAction2::fMax

private

Definition at line 68 of file PrimaryGeneratorAction2.hh.

nPoints

G4int PrimaryGeneratorAction2::nPoints

private

Definition at line 63 of file PrimaryGeneratorAction2.hh.

particleGun

G4ParticleGun* PrimaryGeneratorAction2::particleGun

private

Definition at line 61 of file PrimaryGeneratorAction2.hh.

x

std::vector<G4double> PrimaryGeneratorAction2::x

private

Definition at line 64 of file PrimaryGeneratorAction2.hh.

---

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

• PrimaryGeneratorAction2.hh
• PrimaryGeneratorAction2.cc