#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 37 of file G4HadPhaseSpaceGenbod.hh.

Constructor & Destructor Documentation

◆ G4HadPhaseSpaceGenbod()

G4HadPhaseSpaceGenbod::G4HadPhaseSpaceGenbod ( G4int verbose = 0 )

Definition at line 51 of file G4HadPhaseSpaceGenbod.cc.

52 : G4VHadPhaseSpaceAlgorithm("G4HadPhaseSpaceGenbod",verbose),

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

G4VHadPhaseSpaceAlgorithm::G4VHadPhaseSpaceAlgorithm

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

Definition G4VHadPhaseSpaceAlgorithm.hh:45

◆ ~G4HadPhaseSpaceGenbod()

virtual G4HadPhaseSpaceGenbod::~G4HadPhaseSpaceGenbod ( )

inlinevirtual

Definition at line 40 of file G4HadPhaseSpaceGenbod.hh.

40{;}

Member Function Documentation

◆ AcceptEvent()

G4bool G4HadPhaseSpaceGenbod::AcceptEvent ( ) const

protected

Definition at line 170 of file G4HadPhaseSpaceGenbod.cc.

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

Referenced by GenerateMultiBody().

◆ AccumulateFinalState()

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

protected

Definition at line 196 of file G4HadPhaseSpaceGenbod.cc.

                                                     {
  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;
  }
}

Referenced by GenerateMomenta().

◆ ComputeWeight()

G4double G4HadPhaseSpaceGenbod::ComputeWeight ( ) const

protected

Definition at line 163 of file G4HadPhaseSpaceGenbod.cc.

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

Referenced by AcceptEvent().

◆ ComputeWeightScale()

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

protected

Definition at line 148 of file G4HadPhaseSpaceGenbod.cc.

                                                                           {
  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;
}

Referenced by Initialize().

◆ FillEnergySteps()

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

protected

Definition at line 123 of file G4HadPhaseSpaceGenbod.cc.

                                                          {
  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);
  }
}

Referenced by GenerateMultiBody().

◆ FillRandomBuffer()

void G4HadPhaseSpaceGenbod::FillRandomBuffer ( )

protected

Definition at line 110 of file G4HadPhaseSpaceGenbod.cc.

                                             {
  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);
}

Referenced by GenerateMultiBody().

◆ GenerateMomenta()

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

protected

Definition at line 180 of file G4HadPhaseSpaceGenbod.cc.

                                                {
  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;
  }
}

Referenced by GenerateMultiBody().

◆ GenerateMultiBody()

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

protectedvirtual

Implements G4VHadDecayAlgorithm.

Definition at line 58 of file G4HadPhaseSpaceGenbod.cc.

                                                  {
  if (GetVerboseLevel()) G4cout << GetName() << "::GenerateMultiBody" << G4endl;
 
  finalState.clear();
 
  Initialize(initialMass, masses);
 
  const G4int maxNumberOfLoops = 10000;
  nTrials = 0;
  do {              // Apply accept/reject to get distribution
    ++nTrials;
    FillRandomBuffer();
    FillEnergySteps(initialMass, masses);
  } while ( (!AcceptEvent()) && nTrials < maxNumberOfLoops );  /* Loop checking, 02.11.2015, A.Ribon */
  if ( nTrials >= maxNumberOfLoops ) {
    G4ExceptionDescription ed;
    ed << " Failed sampling after maxNumberOfLoops attempts : forced exit" << G4endl;
    G4Exception( " G4HadPhaseSpaceGenbod::GenerateMultiBody ", "HAD_GENBOD_001", FatalException, ed );
  }
  GenerateMomenta(masses, finalState);
}

◆ Initialize()

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

protected

Definition at line 83 of file G4HadPhaseSpaceGenbod.cc.

                                                                    {
  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);
}

Referenced by GenerateMultiBody().

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

Public Member Functions

Protected Member Functions

Detailed Description

Constructor & Destructor Documentation

◆ G4HadPhaseSpaceGenbod()

◆ ~G4HadPhaseSpaceGenbod()

Member Function Documentation

◆ AcceptEvent()

◆ AccumulateFinalState()

◆ ComputeWeight()

◆ ComputeWeightScale()

◆ FillEnergySteps()

◆ FillRandomBuffer()

◆ GenerateMomenta()

◆ GenerateMultiBody()

◆ Initialize()