#include <G4HadPhaseSpaceGenbod.hh>

Inheritance diagram for G4HadPhaseSpaceGenbod:

Public Member Functions
	G4HadPhaseSpaceGenbod (G4int verbose=0)

virtual	~G4HadPhaseSpaceGenbod ()

Public Member Functions inherited from G4VHadPhaseSpaceAlgorithm
	G4VHadPhaseSpaceAlgorithm (const G4String &algName, G4int verbose=0)

virtual	~G4VHadPhaseSpaceAlgorithm ()

Public Member Functions inherited from G4VHadDecayAlgorithm
	G4VHadDecayAlgorithm (const G4String &algName, G4int verbose=0)

virtual	~G4VHadDecayAlgorithm ()

void	Generate (G4double initialMass, const std::vector< G4double > &masses, std::vector< G4LorentzVector > &finalState)

virtual void	SetVerboseLevel (G4int verbose)

G4int	GetVerboseLevel () const

const G4String &	GetName () const

Protected Member Functions
virtual void	GenerateMultiBody (G4double initialMass, const std::vector< G4double > &masses, std::vector< G4LorentzVector > &finalState)

void	Initialize (G4double initialMass, const std::vector< G4double > &masses)

void	FillRandomBuffer ()

void	ComputeWeightScale (const std::vector< G4double > &masses)

void	FillEnergySteps (G4double initialMass, const std::vector< G4double > &masses)

void	GenerateMomenta (const std::vector< G4double > &masses, std::vector< G4LorentzVector > &finalState)

void	AccumulateFinalState (size_t i, const std::vector< G4double > &masses, std::vector< G4LorentzVector > &finalState)

G4bool	AcceptEvent () const

G4double	ComputeWeight () const

Protected Member Functions inherited from G4VHadPhaseSpaceAlgorithm
virtual void	GenerateTwoBody (G4double initialMass, const std::vector< G4double > &masses, std::vector< G4LorentzVector > &finalState)

G4ThreeVector	UniformVector (G4double mag=1.) const

Protected Member Functions inherited from G4VHadDecayAlgorithm
virtual G4bool	IsDecayAllowed (G4double initialMass, const std::vector< G4double > &masses) const

G4double	TwoBodyMomentum (G4double M0, G4double M1, G4double M2) const

G4double	UniformTheta () const

G4double	UniformPhi () const

void	PrintVector (const std::vector< G4double > &v, const G4String &name, std::ostream &os) const

Detailed Description

Definition at line 38 of file G4HadPhaseSpaceGenbod.hh.

Constructor & Destructor Documentation

G4HadPhaseSpaceGenbod::G4HadPhaseSpaceGenbod ( G4int verbose = 0 )

Definition at line 52 of file G4HadPhaseSpaceGenbod.cc.

53 : G4VHadPhaseSpaceAlgorithm("G4HadPhaseSpaceGenbod",verbose),

54 nFinal(0), totalMass(0.), massExcess(0.), weightMax(0.), nTrials(0) {;}

G4VHadPhaseSpaceAlgorithm::G4VHadPhaseSpaceAlgorithm

G4VHadPhaseSpaceAlgorithm(const G4String &algName, G4int verbose=0)

Definition: G4VHadPhaseSpaceAlgorithm.hh:46

virtual G4HadPhaseSpaceGenbod::~G4HadPhaseSpaceGenbod ( )

inlinevirtual

Definition at line 41 of file G4HadPhaseSpaceGenbod.hh.

41 {;}

Member Function Documentation

G4bool G4HadPhaseSpaceGenbod::AcceptEvent ( ) const

protected

Definition at line 166 of file G4HadPhaseSpaceGenbod.cc.

References ComputeWeight(), G4cout, G4endl, G4UniformRand, G4VHadDecayAlgorithm::GetName(), and G4VHadDecayAlgorithm::GetVerboseLevel().

Referenced by GenerateMultiBody().

                                                 {
   if (GetVerboseLevel()>1) 
     G4cout << GetName() << "::AcceptEvent? " << nTrials << G4endl;
 
   return (G4UniformRand() <= ComputeWeight());
 }

void G4HadPhaseSpaceGenbod::AccumulateFinalState	(	size_t	i,
		const std::vector< G4double > &	masses,
		std::vector< G4LorentzVector > &	finalState
	)

protected

Definition at line 193 of file G4HadPhaseSpaceGenbod.cc.

References G4cout, G4endl, G4UniformRand, G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), and python.hepunit::twopi.

Referenced by GenerateMomenta().

                                                      {
   if (GetVerboseLevel()>2)
     G4cout << GetName() << "::AccumulateFinalState " << i << G4endl;
 
   if (i==0) {           // First final state particle left alone
     finalState[i].setVectM(G4ThreeVector(0.,pd[i],0.),masses[i]);
     return;
   }
   
   finalState[i].setVectM(G4ThreeVector(0.,-pd[i-1],0.),masses[i]);
   G4double phi = G4UniformRand() * twopi;
   G4double theta = std::acos(2.*G4UniformRand() - 1.);
 
   if (GetVerboseLevel() > 2) {
     G4cout << " initialized Py " << -pd[i-1] << " phi " << phi
        << " theta " << theta << G4endl;
   }
 
   G4double esys=0.,beta=0.,gamma=1.;
   if (i < nFinal-1) {           // Do not boost final particle
     esys = std::sqrt(pd[i]*pd[i]+meff[i]*meff[i]);
     beta = pd[i] / esys;
     gamma = esys / meff[i];
 
     if (GetVerboseLevel()>2)
       G4cout << " esys " << esys << " beta " << beta << " gamma " << gamma
          << G4endl;
   }
 
   for (size_t j=0; j<=i; j++) {     // Accumulate rotations
     finalState[j].rotateZ(theta).rotateY(phi);
     finalState[j].setY(gamma*(finalState[j].y() + beta*finalState[j].e()));
     if (GetVerboseLevel()>2) G4cout << " j " << j << " " << finalState[j] << G4endl;
   }
 }

G4double G4HadPhaseSpaceGenbod::ComputeWeight ( ) const

protected

Definition at line 159 of file G4HadPhaseSpaceGenbod.cc.

References G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), and G4VHadDecayAlgorithm::GetVerboseLevel().

Referenced by AcceptEvent().

                                                     {
   if (GetVerboseLevel()>1) G4cout << GetName() << "::ComputeWeight" << G4endl;
 
   return (std::accumulate(pd.begin(), pd.end(), 1./weightMax,
               std::multiplies<G4double>()));
 }

void G4HadPhaseSpaceGenbod::ComputeWeightScale ( const std::vector< G4double > & masses )

protected

Definition at line 144 of file G4HadPhaseSpaceGenbod.cc.

References G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), and G4VHadDecayAlgorithm::TwoBodyMomentum().

Referenced by Initialize().

                                                                            {
   if (GetVerboseLevel()>1) 
     G4cout << GetName() << "::ComputeWeightScale" << G4endl;
 
   weightMax = 1.;
   for (size_t i=1; i<nFinal; i++) {
     weightMax *= TwoBodyMomentum(massExcess+msum[i], msum[i-1], masses[i]);
   }
 
   if (GetVerboseLevel()>2) G4cout << " weightMax = " << weightMax << G4endl;
 }

void G4HadPhaseSpaceGenbod::FillEnergySteps	(	G4double	initialMass,
		const std::vector< G4double > &	masses
	)

protected

Definition at line 119 of file G4HadPhaseSpaceGenbod.cc.

References G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), G4VHadDecayAlgorithm::PrintVector(), and G4VHadDecayAlgorithm::TwoBodyMomentum().

Referenced by GenerateMultiBody().

                                                           {
   if (GetVerboseLevel()>1) G4cout << GetName() << "::FillEnergySteps" << G4endl;
 
   meff.clear();
   pd.clear();
 
   meff.push_back(masses[0]);
   for (size_t i=1; i<nFinal-1; i++) {
     meff.push_back(rndm[i-1]*massExcess + msum[i]);
     pd.push_back(TwoBodyMomentum(meff[i], meff[i-1], masses[i]));
   }
   meff.push_back(initialMass);
   pd.push_back(TwoBodyMomentum(meff[nFinal-1], meff[nFinal-2], masses[nFinal-1]));
 
   if (GetVerboseLevel()>2) {
     PrintVector(meff,"meff",G4cout);
     PrintVector(pd,"pd",G4cout);
   }
 }

void G4HadPhaseSpaceGenbod::FillRandomBuffer ( )

protected

Definition at line 106 of file G4HadPhaseSpaceGenbod.cc.

References G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), G4VHadDecayAlgorithm::PrintVector(), and sort().

Referenced by GenerateMultiBody().

                                              {
   if (GetVerboseLevel()>1) G4cout << GetName() << "::FillRandomBuffer" << G4endl;
 
   rndm.resize(nFinal-2,0.); // Final states generated in sorted order
   std::generate(rndm.begin(), rndm.end(), uniformRand);
   std::sort(rndm.begin(), rndm.end());
   if (GetVerboseLevel()>2) PrintVector(rndm, "rndm", G4cout);
 }

void G4HadPhaseSpaceGenbod::GenerateMomenta	(	const std::vector< G4double > &	masses,
		std::vector< G4LorentzVector > &	finalState
	)

protected

Definition at line 177 of file G4HadPhaseSpaceGenbod.cc.

References AccumulateFinalState(), G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), and G4VHadDecayAlgorithm::GetVerboseLevel().

Referenced by GenerateMultiBody().

                                                 {
   if (GetVerboseLevel()>1) G4cout << GetName() << "::GenerateMomenta" << G4endl;
 
   finalState.resize(nFinal);    // Preallocate vectors for convenience below
 
   for (size_t i=0; i<nFinal; i++) {
     AccumulateFinalState(i, masses, finalState);
     if (GetVerboseLevel()>2)
       G4cout << " finalState[" << i << "] " << finalState[i] << G4endl;
   }
 }

void G4HadPhaseSpaceGenbod::GenerateMultiBody	(	G4double	initialMass,
		const std::vector< G4double > &	masses,
		std::vector< G4LorentzVector > &	finalState
	)

protectedvirtual

Implements G4VHadDecayAlgorithm.

Definition at line 60 of file G4HadPhaseSpaceGenbod.cc.

References AcceptEvent(), FillEnergySteps(), FillRandomBuffer(), G4cout, G4endl, GenerateMomenta(), G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), and Initialize().

                                                   {
   if (GetVerboseLevel()) G4cout << GetName() << "::GenerateMultiBody" << G4endl;
 
   finalState.clear();
 
   Initialize(initialMass, masses);
 
   nTrials = 0;
   do {              // Apply accept/reject to get distribution
     ++nTrials;
     FillRandomBuffer();
     FillEnergySteps(initialMass, masses);
   } while (!AcceptEvent()); // FIXME:  Do we need a limit on nTrials?
 
   GenerateMomenta(masses, finalState);
 }

void G4HadPhaseSpaceGenbod::Initialize	(	G4double	initialMass,
		const std::vector< G4double > &	masses
	)

protected

Definition at line 80 of file G4HadPhaseSpaceGenbod.cc.

References ComputeWeightScale(), G4cout, G4endl, G4VHadDecayAlgorithm::GetName(), G4VHadDecayAlgorithm::GetVerboseLevel(), and G4VHadDecayAlgorithm::PrintVector().

Referenced by GenerateMultiBody().

                                                                     {
   if (GetVerboseLevel()>1) G4cout << GetName() << "::Initialize" << G4endl;
 
   nFinal = masses.size();
   msum.resize(nFinal, 0.);      // Initialize buffers for filling
   msq.resize(nFinal, 0.);
 
   std::partial_sum(masses.begin(), masses.end(), msum.begin());
   std::transform(masses.begin(), masses.end(), masses.begin(), msq.begin(),
          std::multiplies<G4double>());
   totalMass  = msum.back();
   massExcess = initialMass - totalMass;
 
   if (GetVerboseLevel()>2) {
     PrintVector(msum, "msum", G4cout);
     PrintVector(msq, "msq", G4cout);
     G4cout << " totalMass " << totalMass << " massExcess " << massExcess
        << G4endl;
   }
 
   ComputeWeightScale(masses);
 }

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

Public Member Functions

Protected Member Functions

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation