Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CexmcReimplementedGenbod Class Reference

#include <CexmcReimplementedGenbod.hh>

Inheritance diagram for CexmcReimplementedGenbod:
Collaboration diagram for CexmcReimplementedGenbod:

Public Member Functions

 CexmcReimplementedGenbod ()
 
G4double Generate (void)
 
- Public Member Functions inherited from CexmcPhaseSpaceGenerator
 CexmcPhaseSpaceGenerator ()
 
virtual ~CexmcPhaseSpaceGenerator ()
 
virtual G4bool CheckKinematics (void)
 
void SetParticles (const CexmcPhaseSpaceInVector &inVec_, const CexmcPhaseSpaceOutVector &outVec_)
 
void SetFermiEnergyDependence (G4bool on=true)
 

Additional Inherited Members

- Protected Attributes inherited from CexmcPhaseSpaceGenerator
CexmcPhaseSpaceInVector inVec
 
CexmcPhaseSpaceOutVector outVec
 
G4bool fermiEnergyDepIsOn
 
G4double totalEnergy
 
G4double totalMass
 

Detailed Description

Definition at line 51 of file CexmcReimplementedGenbod.hh.

Constructor & Destructor Documentation

CexmcReimplementedGenbod::CexmcReimplementedGenbod ( )

Definition at line 81 of file CexmcReimplementedGenbod.cc.

81  : maxWeight( 0. ),
82  nmbOfOutputParticles( 0 )
83 {
84 }

Member Function Documentation

G4double CexmcReimplementedGenbod::Generate ( void  )
virtual

Implements CexmcPhaseSpaceGenerator.

Definition at line 87 of file CexmcReimplementedGenbod.cc.

88 {
89  // Generate a random final state.
90  // The function returns the weigth of the current event.
91  // Note that Momentum, Energy units are Gev/C, GeV
92 
93  G4double te_minus_tm( totalEnergy - totalMass );
94  G4double rno[ maxParticles ];
95  rno[ 0 ] = 0;
96 
97  if ( nmbOfOutputParticles > 2 )
98  {
99  for ( G4int i( 1 ); i < nmbOfOutputParticles - 1; ++i )
100  {
101  rno[ i ] = G4UniformRand();
102  }
103  qsort( rno + 1, nmbOfOutputParticles - 2, sizeof( G4double ),
104  DoubleMax );
105  }
106  rno[ nmbOfOutputParticles - 1 ] = 1;
107 
108  G4double invMas[ maxParticles ];
109  G4double sum( 0 );
110 
111  for ( int i( 0 ); i < nmbOfOutputParticles; ++i )
112  {
113  sum += outVec[ i ].mass / GeV;
114  invMas[ i ] = rno[ i ] * te_minus_tm / GeV + sum;
115  }
116 
117  //
118  //-----> compute the weight of the current event
119  //
120  G4double wt( maxWeight );
121  G4double pd[ maxParticles ];
122 
123  for ( int i( 0 ); i < nmbOfOutputParticles - 1; ++i )
124  {
125  pd[ i ] = PDK( invMas[ i + 1 ], invMas[ i ],
126  outVec[ i + 1 ].mass / GeV );
127  wt *= pd[ i ];
128  }
129 
130  //
131  //-----> complete specification of event (Raubold-Lynch method)
132  //
133  outVec[ 0 ].lVec->setPx( 0. );
134  outVec[ 0 ].lVec->setPy( pd[ 0 ] );
135  outVec[ 0 ].lVec->setPz( 0. );
136  outVec[ 0 ].lVec->setE( std::sqrt( pd[ 0 ] * pd[ 0 ] +
137  outVec[ 0 ].mass / GeV *
138  outVec[ 0 ].mass / GeV ) );
139 
140  G4int i( 1 );
141 
142  while ( true )
143  {
144  outVec[ i ].lVec->setPx( 0. );
145  outVec[ i ].lVec->setPy( -pd[ i - 1 ] );
146  outVec[ i ].lVec->setPz( 0. );
147  outVec[ i ].lVec->setE( std::sqrt( pd[ i - 1 ] * pd[ i - 1 ] +
148  outVec[ i ].mass / GeV *
149  outVec[ i ].mass / GeV ) );
150 
151  G4double cZ( 2 * G4UniformRand() - 1 );
152  G4double sZ( std::sqrt( 1 - cZ * cZ ) );
153  G4double angY( 2 * pi * G4UniformRand() );
154  G4double cY( std::cos( angY ) );
155  G4double sY( std::sin( angY ) );
156 
157  for ( int j( 0 ); j <= i; ++j )
158  {
159  G4LorentzVector * v( outVec[ j ].lVec );
160  G4double x( v->px() );
161  G4double y( v->py() );
162  v->setPx( cZ * x - sZ * y );
163  v->setPy( sZ * x + cZ * y ); // rotation around Z
164  x = v->px();
165  G4double z( v->pz() );
166  v->setPx( cY * x - sY * z );
167  v->setPz( sY * x + cY * z ); // rotation around Y
168  }
169 
170  if ( i == nmbOfOutputParticles - 1 )
171  break;
172 
173  G4double beta( pd[ i ] / std::sqrt( pd[ i ] * pd[ i ] +
174  invMas[ i ] * invMas[ i ] ) );
175  for ( int j( 0 ); j <= i; ++j )
176  outVec[ j ].lVec->boost( 0, beta, 0 );
177 
178  ++i;
179  }
180 
181  for ( int i( 0 ); i < nmbOfOutputParticles; ++i )
182  *outVec[ i ].lVec *= GeV;
183 
184  //
185  //---> return the weigth of event
186  //
187  return wt;
188 }
tuple x
Definition: test.py:50
CexmcPhaseSpaceOutVector outVec
int G4int
Definition: G4Types.hh:78
#define G4UniformRand()
Definition: Randomize.hh:97
tuple v
Definition: test.py:18
static constexpr double GeV
Definition: G4SIunits.hh:217
tuple z
Definition: test.py:28
static constexpr double pi
Definition: G4SIunits.hh:75
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:


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