Geant4-11
Data Structures | Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes
G4ParticleHPFSFissionFS Class Reference

#include <G4ParticleHPFSFissionFS.hh>

Inheritance diagram for G4ParticleHPFSFissionFS:
G4ParticleHPFinalState

Data Structures

struct  toBeCached
 

Public Member Functions

G4DynamicParticleVectorApplyYourself (G4int Prompt, G4int delayed, G4double *decayconst)
 
 G4ParticleHPFSFissionFS ()
 
G4double GetA ()
 
G4ParticleHPFissionEReleaseGetEnergyRelease ()
 
G4int GetM ()
 
G4double GetMass ()
 
G4double GetN ()
 
G4DynamicParticleVectorGetPhotons ()
 
virtual G4ParticleHPVectorGetXsec ()
 
virtual G4double GetXsec (G4double)
 
G4double GetZ ()
 
G4bool HasAnyData ()
 
G4bool HasFSData ()
 
G4bool HasXsec ()
 
void Init (G4double A, G4double Z, G4int M, G4String &dirName, G4String &aFSType, G4ParticleDefinition *)
 
void Init (G4double A, G4double Z, G4String &dirName, G4String &aFSType, G4ParticleDefinition *projectile)
 
G4ParticleHPFinalStateNew ()
 
void SampleNeutronMult (G4int &all, G4int &Prompt, G4int &delayed, G4double energy, G4int off)
 
void SetA_Z (G4double anA, G4double aZ, G4int aM=0)
 
void SetAZMs (G4double anA, G4double aZ, G4int aM, G4ParticleHPDataUsed used)
 
void SetNeutronRP (const G4ReactionProduct &aNeutron)
 
void SetProjectile (G4ParticleDefinition *projectile)
 
void SetTarget (const G4ReactionProduct &aTarget)
 
 ~G4ParticleHPFSFissionFS ()
 

Protected Member Functions

void adjust_final_state (G4LorentzVector)
 

Protected Attributes

G4bool hasAnyData
 
G4bool hasFSData
 
G4bool hasXsec
 
G4int secID
 
G4double theBaseA
 
G4int theBaseM
 
G4double theBaseZ
 
G4int theNDLDataA
 
G4int theNDLDataM
 
G4int theNDLDataZ
 
G4ParticleDefinitiontheProjectile
 
G4Cache< G4HadFinalState * > theResult
 

Private Member Functions

G4HadFinalStateApplyYourself (const G4HadProjectile &)
 

Private Attributes

G4Cache< toBeCachedfCache
 
G4ParticleHPEnergyDistribution theDelayedNeutronEnDis
 
G4ParticleHPFissionERelease theEnergyRelease
 
G4ParticleHPParticleYield theFinalStateNeutrons
 
G4ParticleHPPhotonDist theFinalStatePhotons
 
G4ParticleHPNames theNames
 
G4ParticleHPAngular theNeutronAngularDis
 
G4ParticleHPEnergyDistribution thePromptNeutronEnDis
 

Detailed Description

Definition at line 45 of file G4ParticleHPFSFissionFS.hh.

Constructor & Destructor Documentation

◆ G4ParticleHPFSFissionFS()

G4ParticleHPFSFissionFS::G4ParticleHPFSFissionFS ( )
inline

Definition at line 56 of file G4ParticleHPFSFissionFS.hh.

References G4ParticleHPFinalState::hasXsec.

Referenced by New().

◆ ~G4ParticleHPFSFissionFS()

G4ParticleHPFSFissionFS::~G4ParticleHPFSFissionFS ( )
inline

Definition at line 57 of file G4ParticleHPFSFissionFS.hh.

57{}

Member Function Documentation

◆ adjust_final_state()

void G4ParticleHPFinalState::adjust_final_state ( G4LorentzVector  init_4p_lab)
protectedinherited

Definition at line 47 of file G4ParticleHPFinalState.cc.

48{
49
50 G4double minimum_energy = 1*keV;
51
52 if ( G4ParticleHPManager::GetInstance()->GetDoNotAdjustFinalState() ) return;
53
54 G4int nSecondaries = theResult.Get()->GetNumberOfSecondaries();
55
56 G4int sum_Z = 0;
57 G4int sum_A = 0;
58 G4int max_SecZ = 0;
59 G4int max_SecA = 0;
60 G4int imaxA = -1;
61 for ( int i = 0 ; i < nSecondaries ; i++ )
62 {
63 //G4cout << "G4ParticleHPFinalState::adjust_final_state theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetParticleName() = " << theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetParticleName() << G4endl;
65 max_SecZ = std::max ( max_SecZ , theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetAtomicNumber() );
67 max_SecA = std::max ( max_SecA , theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetAtomicMass() );
68 if ( theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetAtomicMass() == max_SecA ) imaxA = i;
69#ifdef G4PHPDEBUG
70 if( std::getenv("G4ParticleHPDebug")) G4cout << "G4ParticleHPFinalState::adjust_final_stat SECO " << i << " " <<theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetParticleName() << G4endl;
71#endif
72
73 }
74
75 G4ParticleDefinition* resi_pd = 0;
76
77 G4double baseZNew = theBaseZ;
78 G4double baseANew = theBaseA;
80 baseANew ++;
81 } else if( theProjectile == G4Proton::Proton() ) {
82 baseZNew ++;
83 baseANew ++;
84 } else if( theProjectile == G4Deuteron::Deuteron() ) {
85 baseZNew ++;
86 baseANew += 2;
87 } else if( theProjectile == G4Triton::Triton() ) {
88 baseZNew ++;
89 baseANew += 3;
90 } else if( theProjectile == G4He3::He3() ) {
91 baseZNew += 2;
92 baseANew += 3;
93 } else if( theProjectile == G4Alpha::Alpha() ) {
94 baseZNew += 2;
95 baseANew += 4;
96 }
97
98#ifdef G4PHPDEBUG
99 if( std::getenv("G4ParticleHPDebug")) G4cout << "G4ParticleHPFinalState::adjust_final_stat BaseZ " << baseZNew << " BaseA " << baseANew << " sum_Z " << sum_Z << " sum_A " << sum_A << G4endl;
100#endif
101
102 G4bool needOneMoreSec = false;
103 G4ParticleDefinition* oneMoreSec_pd = 0;
104 if ( (int)(baseZNew - sum_Z) == 0 && (int)(baseANew - sum_A) == 0 )
105 {
106 //All secondaries are already created;
107 resi_pd = theResult.Get()->GetSecondary( imaxA )->GetParticle()->GetDefinition();
108 }
109 else
110 {
111 if ( max_SecA >= int(baseANew - sum_A) )
112 {
113 //Most heavy secondary is interpreted as residual
114 resi_pd = theResult.Get()->GetSecondary( imaxA )->GetParticle()->GetDefinition();
115 needOneMoreSec = true;
116 }
117 else
118 {
119 //creation of residual is requierd
120 resi_pd = G4IonTable::GetIonTable()->GetIon ( int(baseZNew - sum_Z) , (int)(baseANew - sum_A) , 0.0 );
121 }
122
123 if ( needOneMoreSec )
124 {
125 if ( int(baseZNew - sum_Z) == 0 && (int)(baseANew - sum_A) > 0 )
126 {
127 //In this case, one neutron is added to secondaries
128 if ( int(baseANew - sum_A) > 1 ) G4cout << "More than one neutron is required for the balance of baryon number!" << G4endl;
129 oneMoreSec_pd = G4Neutron::Neutron();
130 }
131 else
132 {
133#ifdef G4PHPDEBUG
134 if( std::getenv("G4ParticleHPDebug")) G4cout << this << "G4ParticleHPFinalState oneMoreSec_pd Z " << baseZNew << " - " << sum_Z << " A " << baseANew << " - " << sum_A << " projectile " << theProjectile->GetParticleName() << G4endl;
135#endif
136 oneMoreSec_pd = G4IonTable::GetIonTable()->GetIon ( int(baseZNew - sum_Z) , (int)(baseANew - sum_A) , 0.0 );
137 if( !oneMoreSec_pd ) {
138 G4cerr << this << "G4ParticleHPFinalState oneMoreSec_pd Z " << baseZNew << " - " << sum_Z << " A " << baseANew << " - " << sum_A << " projectile " << theProjectile->GetParticleName() << G4endl;
139 G4Exception("G4ParticleHPFinalState:adjust_final_state",
140 "Warning",
142 "No adjustment will be done!");
143 return;
144 }
145 }
146 }
147
148 if ( resi_pd == 0 )
149 {
150 // theNDLDataZ,A has the Z and A of used NDL file
151 G4double ndlZNew = theNDLDataZ;
152 G4double ndlANew = theNDLDataA;
154 ndlANew ++;
155 } else if( theProjectile == G4Proton::Proton() ) {
156 ndlZNew ++;
157 ndlANew ++;
158 } else if( theProjectile == G4Deuteron::Deuteron() ) {
159 ndlZNew ++;
160 ndlANew += 2;
161 } else if( theProjectile == G4Triton::Triton() ) {
162 ndlZNew ++;
163 ndlANew += 3;
164 } else if( theProjectile == G4He3::He3() ) {
165 ndlZNew += 2;
166 ndlANew += 3;
167 } else if( theProjectile == G4Alpha::Alpha() ) {
168 ndlZNew += 2;
169 ndlANew += 4;
170 }
171 // theNDLDataZ,A has the Z and A of used NDL file
172 if ( (int)(ndlZNew - sum_Z) == 0 && (int)(ndlANew - sum_A) == 0 )
173 {
174 G4int dif_Z = ( int ) ( theNDLDataZ - theBaseZ );
175 G4int dif_A = ( int ) ( theNDLDataA - theBaseA );
176 resi_pd = G4IonTable::GetIonTable()->GetIon ( max_SecZ - dif_Z , max_SecA - dif_A , 0.0 );
177 if( !resi_pd ) {
178 G4cerr << "G4ParticleHPFinalState resi_pd Z " << max_SecZ << " - " << dif_Z << " A " << max_SecA << " - " << dif_A << " projectile " << theProjectile->GetParticleName() << G4endl;
179 G4Exception("G4ParticleHPFinalState:adjust_final_state",
180 "Warning",
182 "No adjustment will be done!");
183 return;
184 }
185
186 for ( int i = 0 ; i < nSecondaries ; i++ )
187 {
188 if ( theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetAtomicNumber() == max_SecZ
189 && theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition()->GetAtomicMass() == max_SecA )
190 {
192 p = p * resi_pd->GetPDGMass()/ G4IonTable::GetIonTable()->GetIon ( max_SecZ , max_SecA , 0.0 )->GetPDGMass();
193 theResult.Get()->GetSecondary( i )->GetParticle()->SetDefinition( resi_pd );
195 }
196 }
197 }
198 }
199 }
200
201
202 G4LorentzVector secs_4p_lab( 0.0 );
203
205 G4double fast = 0;
206 G4double slow = 1;
207 G4int ifast = 0;
208 G4int islow = 0;
209 G4int ires = -1;
210
211 for ( G4int i = 0 ; i < n_sec ; i++ )
212 {
213
214 //G4cout << "HP_DB " << i
215 // << " " << theResult.GetSecondary( i )->GetParticle()->GetDefinition()->GetParticleName()
216 // << " 4p " << theResult.GetSecondary( i )->GetParticle()->Get4Momentum()
217 // << " ke " << theResult.GetSecondary( i )->GetParticle()->Get4Momentum().e() - theResult.GetSecondary( i )->GetParticle()->GetDefinition()->GetPDGMass()
218 // << G4endl;
219
220 secs_4p_lab += theResult.Get()->GetSecondary( i )->GetParticle()->Get4Momentum();
221
222 G4double beta = 0;
224 {
226 }
227 else
228 {
229 beta = 1;
230 }
231
232 if ( theResult.Get()->GetSecondary( i )->GetParticle()->GetDefinition() == resi_pd ) ires = i;
233
234 if ( slow > beta && beta != 0 )
235 {
236 slow = beta;
237 islow = i;
238 }
239
240 if ( fast <= beta )
241 {
242 if ( fast != 1 )
243 {
244 fast = beta;
245 ifast = i;
246 }
247 else
248 {
249// fast is already photon then check E
251 if ( e > theResult.Get()->GetSecondary( ifast )->GetParticle()->Get4Momentum().e() )
252 {
253// among photons, the highest E becomes the fastest
254 ifast = i;
255 }
256 }
257 }
258 }
259
260
261 G4LorentzVector dif_4p = init_4p_lab - secs_4p_lab;
262
263 //G4cout << "HP_DB dif_4p " << init_4p_lab - secs_4p_lab << G4endl;
264 //G4cout << "HP_DB dif_3p mag " << ( dif_4p.v() ).mag() << G4endl;
265 //G4cout << "HP_DB dif_e " << dif_4p.e() - ( dif_4p.v() ).mag()<< G4endl;
266
267 G4LorentzVector p4(0);
268 if ( ires == -1 )
269 {
270// Create and Add Residual Nucleus
271 ires = nSecondaries;
272 nSecondaries += 1;
273
274 G4DynamicParticle* res = new G4DynamicParticle ( resi_pd , dif_4p.v() );
275 theResult.Get()->AddSecondary ( res, secID );
276
277 p4 = res->Get4Momentum();
278 if ( slow > p4.beta() )
279 {
280 slow = p4.beta();
281 islow = ires;
282 }
283 dif_4p = init_4p_lab - ( secs_4p_lab + p4 );
284 }
285
286 if ( needOneMoreSec && oneMoreSec_pd)
287 //
288 // fca: this is not a fix, this is a crash avoidance...
289 // fca: the baryon number is still wrong, most probably because it
290 // fca: should have been decreased, but since we could not create a particle
291 // fca: we just do not add it
292 //
293 {
294 nSecondaries += 1;
295 G4DynamicParticle* one = new G4DynamicParticle ( oneMoreSec_pd , dif_4p.v() );
296 theResult.Get()->AddSecondary ( one, secID );
297 p4 = one->Get4Momentum();
298 if ( slow > p4.beta() )
299 {
300 slow = p4.beta();
301 islow = nSecondaries-1; //Because the first is 0th, so the last becomes "nSecondaries-1"
302 }
303 dif_4p = init_4p_lab - ( secs_4p_lab + p4 );
304 }
305
306 //Which is bigger dif_p or dif_e
307
308 if ( dif_4p.v().mag() < std::abs( dif_4p.e() ) )
309 {
310
311 // Adjust p
312 //if ( dif_4p.v().mag() < 1*MeV )
313 if ( minimum_energy < dif_4p.v().mag() && dif_4p.v().mag() < 1*MeV )
314 {
315
316 nSecondaries += 1;
317 theResult.Get()->AddSecondary ( new G4DynamicParticle ( G4Gamma::Gamma() , dif_4p.v() ), secID );
318
319 }
320 else
321 {
322 //G4cout << "HP_DB Difference in dif_p is too large (>1MeV) or too small(<1keV) to adjust, so that give up tuning" << G4endl;
323 }
324
325 }
326 else
327 {
328
329 // dif_p > dif_e
330 // at first momentum
331 // Move residual momentum
332
333 p4 = theResult.Get()->GetSecondary( ires )->GetParticle()->Get4Momentum();
334 theResult.Get()->GetSecondary( ires )->GetParticle()->SetMomentum( p4.v() + dif_4p.v() );
335 dif_4p = init_4p_lab - ( secs_4p_lab - p4 + theResult.Get()->GetSecondary( ires )->GetParticle()->Get4Momentum() );
336
337 //G4cout << "HP_DB new residual kinetic energy " << theResult.GetSecondary( ires )->GetParticle()->GetKineticEnergy() << G4endl;
338
339 }
340
341 G4double dif_e = dif_4p.e() - ( dif_4p.v() ).mag();
342 //G4cout << "HP_DB dif_e " << dif_e << G4endl;
343
344 if ( dif_e > 0 )
345 {
346
347// create 2 gamma
348
349 nSecondaries += 2;
350 G4double e1 = ( dif_4p.e() -dif_4p.v().mag() ) / 2;
351
352 if ( minimum_energy < e1 )
353 {
354 G4double costh = 2.*G4UniformRand()-1.;
356 G4ThreeVector dir( std::sin(std::acos(costh))*std::cos(phi),
357 std::sin(std::acos(costh))*std::sin(phi),
358 costh);
361 }
362 else
363 {
364 //G4cout << "HP_DB Difference is too small(<1keV) to adjust, so that neglect it" << G4endl;
365 }
366
367 }
368 else //dif_e < 0
369 {
370
371// At first reduce KE of the fastest secondary;
374 G4ThreeVector dir = ( theResult.Get()->GetSecondary( ifast )->GetParticle()->GetMomentum() ).unit();
375
376 //G4cout << "HP_DB ifast " << ifast << " ke0 " << ke0 << G4endl;
377
378 if ( ke0 + dif_e > 0 )
379 {
380 theResult.Get()->GetSecondary( ifast )->GetParticle()->SetKineticEnergy( ke0 + dif_e );
381 G4ThreeVector dp = p0 - theResult.Get()->GetSecondary( ifast )->GetParticle()->GetMomentum();
382
384 //theResult.GetSecondary( islow )->GetParticle()->SetMomentum( p - dif_e*dir );
385 theResult.Get()->GetSecondary( islow )->GetParticle()->SetMomentum( p + dp );
386 }
387 else
388 {
389 //G4cout << "HP_DB Difference in dif_e too large ( <0MeV ) to adjust, so that give up tuning" << G4endl;
390 }
391
392 }
393
394}
static const G4double e1[44]
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
static constexpr double twopi
Definition: G4SIunits.hh:56
static constexpr double keV
Definition: G4SIunits.hh:202
static constexpr double MeV
Definition: G4SIunits.hh:200
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
G4GLOB_DLL std::ostream G4cerr
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
#define G4UniformRand()
Definition: Randomize.hh:52
double mag() const
Hep3Vector v() const
static G4Alpha * Alpha()
Definition: G4Alpha.cc:88
value_type & Get() const
Definition: G4Cache.hh:315
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:93
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4ParticleDefinition * GetDefinition() const
G4LorentzVector Get4Momentum() const
G4double GetKineticEnergy() const
void SetMomentum(const G4ThreeVector &momentum)
G4ThreeVector GetMomentum() const
void SetKineticEnergy(G4double aEnergy)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
std::size_t GetNumberOfSecondaries() const
G4HadSecondary * GetSecondary(size_t i)
G4DynamicParticle * GetParticle()
static G4He3 * He3()
Definition: G4He3.cc:93
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:522
static G4IonTable * GetIonTable()
Definition: G4IonTable.cc:170
static G4Neutron * Neutron()
Definition: G4Neutron.cc:103
G4int GetAtomicNumber() const
G4int GetAtomicMass() const
const G4String & GetParticleName() const
G4ParticleDefinition * theProjectile
G4Cache< G4HadFinalState * > theResult
static G4ParticleHPManager * GetInstance()
static G4Proton * Proton()
Definition: G4Proton.cc:92
static G4Triton * Triton()
Definition: G4Triton.cc:93
T max(const T t1, const T t2)
brief Return the largest of the two arguments

References G4HadFinalState::AddSecondary(), G4Alpha::Alpha(), CLHEP::HepLorentzVector::beta(), anonymous_namespace{G4PionRadiativeDecayChannel.cc}::beta, G4Deuteron::Deuteron(), CLHEP::HepLorentzVector::e(), e1, G4cerr, G4cout, G4endl, G4Exception(), G4UniformRand, G4Gamma::Gamma(), G4Cache< VALTYPE >::Get(), G4DynamicParticle::Get4Momentum(), G4ParticleDefinition::GetAtomicMass(), G4ParticleDefinition::GetAtomicNumber(), G4DynamicParticle::GetDefinition(), G4ParticleHPManager::GetInstance(), G4IonTable::GetIon(), G4IonTable::GetIonTable(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentum(), G4HadFinalState::GetNumberOfSecondaries(), G4HadSecondary::GetParticle(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGMass(), G4HadFinalState::GetSecondary(), G4He3::He3(), JustWarning, keV, CLHEP::Hep3Vector::mag(), G4INCL::Math::max(), MeV, G4Neutron::Neutron(), G4Proton::Proton(), G4ParticleHPFinalState::secID, G4DynamicParticle::SetDefinition(), G4DynamicParticle::SetKineticEnergy(), G4DynamicParticle::SetMomentum(), G4ParticleHPFinalState::theBaseA, G4ParticleHPFinalState::theBaseZ, G4ParticleHPFinalState::theNDLDataA, G4ParticleHPFinalState::theNDLDataZ, G4ParticleHPFinalState::theProjectile, G4ParticleHPFinalState::theResult, G4Triton::Triton(), twopi, and CLHEP::HepLorentzVector::v().

Referenced by G4ParticleHPInelasticBaseFS::BaseApply(), and G4ParticleHPInelasticCompFS::CompositeApply().

◆ ApplyYourself() [1/2]

G4HadFinalState * G4ParticleHPFSFissionFS::ApplyYourself ( const G4HadProjectile )
inlineprivatevirtual

Reimplemented from G4ParticleHPFinalState.

Definition at line 103 of file G4ParticleHPFSFissionFS.hh.

103{ return 0; }

◆ ApplyYourself() [2/2]

G4DynamicParticleVector * G4ParticleHPFSFissionFS::ApplyYourself ( G4int  Prompt,
G4int  delayed,
G4double decayconst 
)

Definition at line 103 of file G4ParticleHPFSFissionFS.cc.

105 {
106 G4int i;
108 G4ReactionProduct boosted;
109 boosted.Lorentz( *(fCache.Get().theNeutronRP) , *(fCache.Get().theTarget) );
110 G4double eKinetic = boosted.GetKineticEnergy();
111
112// Build neutrons
113 G4ReactionProduct * theNeutrons = new G4ReactionProduct[nPrompt+nDelayed];
114 for(i=0; i<nPrompt+nDelayed; i++)
115 {
116 theNeutrons[i].SetDefinition(G4Neutron::Neutron());
117 }
118
119// sample energies
120 G4int it, dummy;
121 G4double tempE;
122 for(i=0; i<nPrompt; i++)
123 {
124 tempE = thePromptNeutronEnDis.Sample(eKinetic, dummy); // energy distribution (file5) always in lab
125 theNeutrons[i].SetKineticEnergy(tempE);
126 }
127 for(i=nPrompt; i<nPrompt+nDelayed; i++)
128 {
129 theNeutrons[i].SetKineticEnergy(theDelayedNeutronEnDis.Sample(eKinetic, it)); // dito
130 if(it==0) theNeutrons[i].SetKineticEnergy(thePromptNeutronEnDis.Sample(eKinetic, dummy));
131 theDecayConst[i-nPrompt] = theFinalStateNeutrons.GetDecayConstant(it); // this is returned
132 }
133
134// sample neutron angular distribution
135 for(i=0; i<nPrompt+nDelayed; i++)
136 {
137 theNeutronAngularDis.SampleAndUpdate(theNeutrons[i]); // angular comes back in lab automatically
138 }
139
140// already in lab. Add neutrons to dynamic particle vector
141 for(i=0; i<nPrompt+nDelayed; i++)
142 {
144 dp->SetDefinition(theNeutrons[i].GetDefinition());
145 dp->SetMomentum(theNeutrons[i].GetMomentum());
146 aResult->push_back(dp);
147 }
148 delete [] theNeutrons;
149// return the result
150 return aResult;
151 }
std::vector< G4DynamicParticle * > G4DynamicParticleVector
void SampleAndUpdate(G4ReactionProduct &anIncidentParticle)
G4double Sample(G4double anEnergy, G4int &it)
G4ParticleHPAngular theNeutronAngularDis
G4ParticleHPParticleYield theFinalStateNeutrons
G4ParticleHPEnergyDistribution thePromptNeutronEnDis
G4ParticleHPEnergyDistribution theDelayedNeutronEnDis
G4double GetKineticEnergy() const
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
void SetKineticEnergy(const G4double en)

References fCache, G4ParticleHPParticleYield::GetDecayConstant(), G4ReactionProduct::GetKineticEnergy(), G4ReactionProduct::Lorentz(), G4Neutron::Neutron(), G4ParticleHPEnergyDistribution::Sample(), G4ParticleHPAngular::SampleAndUpdate(), G4DynamicParticle::SetDefinition(), G4ReactionProduct::SetDefinition(), G4ReactionProduct::SetKineticEnergy(), G4DynamicParticle::SetMomentum(), theDelayedNeutronEnDis, theFinalStateNeutrons, theNeutronAngularDis, and thePromptNeutronEnDis.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

◆ GetA()

G4double G4ParticleHPFinalState::GetA ( void  )
inlineinherited

Definition at line 104 of file G4ParticleHPFinalState.hh.

104{ return theBaseA; }

References G4ParticleHPFinalState::theBaseA.

◆ GetEnergyRelease()

G4ParticleHPFissionERelease * G4ParticleHPFSFissionFS::GetEnergyRelease ( )
inline

Definition at line 96 of file G4ParticleHPFSFissionFS.hh.

97 {
98 return &theEnergyRelease;
99 }
G4ParticleHPFissionERelease theEnergyRelease

References theEnergyRelease.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

◆ GetM()

G4int G4ParticleHPFinalState::GetM ( )
inlineinherited

◆ GetMass()

G4double G4ParticleHPFSFissionFS::GetMass ( )
inline

◆ GetN()

G4double G4ParticleHPFinalState::GetN ( )
inlineinherited

◆ GetPhotons()

G4DynamicParticleVector * G4ParticleHPFSFissionFS::GetPhotons ( )

Definition at line 177 of file G4ParticleHPFSFissionFS.cc.

178{
179// sample photons
181 G4ReactionProduct boosted;
182// the photon distributions are in the Nucleus rest frame.
183 boosted.Lorentz( *(fCache.Get().theNeutronRP) , *(fCache.Get().theTarget) );
184 G4double anEnergy = boosted.GetKineticEnergy();
185 temp = theFinalStatePhotons.GetPhotons(anEnergy);
186 if(temp == 0) { return 0; }
187
188// lorentz transform, and add photons to final state
189 unsigned int i;
191 for(i=0; i<temp->size(); i++)
192 {
193 // back to lab
194 temp->operator[](i)->Lorentz(*(temp->operator[](i)), -1.* (*(fCache.Get().theTarget)) );
196 theOne->SetDefinition(temp->operator[](i)->GetDefinition());
197 theOne->SetMomentum(temp->operator[](i)->GetMomentum());
198 result->push_back(theOne);
199 delete temp->operator[](i);
200 }
201 delete temp;
202 return result;
203}
std::vector< G4ReactionProduct * > G4ReactionProductVector
G4ParticleHPPhotonDist theFinalStatePhotons
G4ReactionProductVector * GetPhotons(G4double anEnergy)

References fCache, G4ReactionProduct::GetKineticEnergy(), G4ParticleHPPhotonDist::GetPhotons(), G4ReactionProduct::Lorentz(), G4DynamicParticle::SetDefinition(), G4DynamicParticle::SetMomentum(), and theFinalStatePhotons.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

◆ GetXsec() [1/2]

virtual G4ParticleHPVector * G4ParticleHPFinalState::GetXsec ( )
inlinevirtualinherited

Reimplemented in G4ParticleHPFissionBaseFS, G4ParticleHPInelasticBaseFS, and G4ParticleHPInelasticCompFS.

Definition at line 99 of file G4ParticleHPFinalState.hh.

99{ return 0; }

◆ GetXsec() [2/2]

virtual G4double G4ParticleHPFinalState::GetXsec ( G4double  )
inlinevirtualinherited

◆ GetZ()

G4double G4ParticleHPFinalState::GetZ ( void  )
inlineinherited

◆ HasAnyData()

G4bool G4ParticleHPFinalState::HasAnyData ( )
inlineinherited

◆ HasFSData()

G4bool G4ParticleHPFinalState::HasFSData ( )
inlineinherited

◆ HasXsec()

G4bool G4ParticleHPFinalState::HasXsec ( )
inlineinherited

Definition at line 94 of file G4ParticleHPFinalState.hh.

94{ return hasXsec; }

References G4ParticleHPFinalState::hasXsec.

Referenced by G4ParticleHPChannel::DumpInfo().

◆ Init() [1/2]

void G4ParticleHPFSFissionFS::Init ( G4double  A,
G4double  Z,
G4int  M,
G4String dirName,
G4String aFSType,
G4ParticleDefinition  
)
virtual

Implements G4ParticleHPFinalState.

Definition at line 45 of file G4ParticleHPFSFissionFS.cc.

46 {
47 G4String tString = "/FS/";
48 G4bool dbool;
49 G4ParticleHPDataUsed aFile = theNames.GetName(static_cast<G4int>(A), static_cast<G4int>(Z), M, dirName, tString, dbool);
50 G4String filename = aFile.GetName();
51 SetAZMs( A, Z, M, aFile );
52 if(!dbool)
53 {
54 hasAnyData = false;
55 hasFSData = false;
56 hasXsec = false;
57 return;
58 }
59 //std::ifstream theData(filename, std::ios::in);
60 std::istringstream theData(std::ios::in);
62
63 // here it comes
64 G4int infoType, dataType;
65 hasFSData = false;
66 while (theData >> infoType) // Loop checking, 11.05.2015, T. Koi
67 {
68 hasFSData = true;
69 theData >> dataType;
70 switch(infoType)
71 {
72 case 1:
73 if(dataType==4) theNeutronAngularDis.Init(theData);
74 if(dataType==5) thePromptNeutronEnDis.Init(theData);
75 if(dataType==12) theFinalStatePhotons.InitMean(theData);
76 if(dataType==14) theFinalStatePhotons.InitAngular(theData);
77 if(dataType==15) theFinalStatePhotons.InitEnergies(theData);
78 break;
79 case 2:
80 if(dataType==1) theFinalStateNeutrons.InitMean(theData);
81 break;
82 case 3:
83 if(dataType==1) theFinalStateNeutrons.InitDelayed(theData);
84 if(dataType==5) theDelayedNeutronEnDis.Init(theData);
85 break;
86 case 4:
87 if(dataType==1) theFinalStateNeutrons.InitPrompt(theData);
88 break;
89 case 5:
90 if(dataType==1) theEnergyRelease.Init(theData);
91 break;
92 default:
93 G4cout << "G4ParticleHPFSFissionFS::Init: unknown data type"<<dataType<<G4endl;
94 throw G4HadronicException(__FILE__, __LINE__, "G4ParticleHPFSFissionFS::Init: unknown data type");
95 break;
96 }
97 }
98 //targetMass = theFinalStateNeutrons.GetTargetMass();
99 //theData.close();
100 }
#define M(row, col)
const G4int Z[17]
const G4double A[17]
void Init(std::istream &aDataFile)
void SetAZMs(G4double anA, G4double aZ, G4int aM, G4ParticleHPDataUsed used)
void Init(std::istream &aDataFile)
void GetDataStream(G4String, std::istringstream &iss)
G4ParticleHPDataUsed GetName(G4int A, G4int Z, G4String base, G4String rest, G4bool &active)
void InitMean(std::istream &aDataFile)
void InitDelayed(std::istream &aDataFile)
void InitPrompt(std::istream &aDataFile)
void InitEnergies(std::istream &aDataFile)
void InitAngular(std::istream &aDataFile)
G4bool InitMean(std::istream &aDataFile)

References A, G4cout, G4endl, G4ParticleHPManager::GetDataStream(), G4ParticleHPManager::GetInstance(), G4ParticleHPDataUsed::GetName(), G4ParticleHPNames::GetName(), G4ParticleHPFinalState::hasAnyData, G4ParticleHPFinalState::hasFSData, G4ParticleHPFinalState::hasXsec, G4ParticleHPAngular::Init(), G4ParticleHPFissionERelease::Init(), G4ParticleHPEnergyDistribution::Init(), G4ParticleHPPhotonDist::InitAngular(), G4ParticleHPParticleYield::InitDelayed(), G4ParticleHPPhotonDist::InitEnergies(), G4ParticleHPParticleYield::InitMean(), G4ParticleHPPhotonDist::InitMean(), G4ParticleHPParticleYield::InitPrompt(), M, G4ParticleHPFinalState::SetAZMs(), theDelayedNeutronEnDis, theEnergyRelease, theFinalStateNeutrons, theFinalStatePhotons, theNames, theNeutronAngularDis, thePromptNeutronEnDis, and Z.

Referenced by G4ParticleHPFissionFS::Init().

◆ Init() [2/2]

void G4ParticleHPFinalState::Init ( G4double  A,
G4double  Z,
G4String dirName,
G4String aFSType,
G4ParticleDefinition projectile 
)
inlineinherited

Definition at line 76 of file G4ParticleHPFinalState.hh.

78 {
79 G4int M = 0;
80 Init ( A, Z, M, dirName, aFSType,const_cast<G4ParticleDefinition*>(projectile));
81 }
void Init(G4double A, G4double Z, G4String &dirName, G4String &aFSType, G4ParticleDefinition *projectile)

References A, G4ParticleHPFinalState::Init(), M, and Z.

Referenced by G4ParticleHPFinalState::Init(), and G4ParticleHPChannel::UpdateData().

◆ New()

G4ParticleHPFinalState * G4ParticleHPFSFissionFS::New ( )
inlinevirtual

◆ SampleNeutronMult()

void G4ParticleHPFSFissionFS::SampleNeutronMult ( G4int all,
G4int Prompt,
G4int delayed,
G4double  energy,
G4int  off 
)

Definition at line 153 of file G4ParticleHPFSFissionFS.cc.

154{
155 G4double promptNeutronMulti = 0;
156 promptNeutronMulti = theFinalStateNeutrons.GetPrompt(eKinetic);
157 G4double delayedNeutronMulti = 0;
158 delayedNeutronMulti = theFinalStateNeutrons.GetDelayed(eKinetic);
159
160 if(delayedNeutronMulti==0&&promptNeutronMulti==0)
161 {
162 Prompt = 0;
163 delayed = 0;
164 G4double totalNeutronMulti = theFinalStateNeutrons.GetMean(eKinetic);
165 all = G4Poisson(totalNeutronMulti-off);
166 all += off;
167 }
168 else
169 {
170 Prompt = G4Poisson(promptNeutronMulti-off);
171 Prompt += off;
172 delayed = G4Poisson(delayedNeutronMulti);
173 all = Prompt+delayed;
174 }
175}
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:50
G4double GetMean(G4double anEnergy)
G4double GetPrompt(G4double anEnergy)
G4double GetDelayed(G4double anEnergy)

References G4Poisson(), G4ParticleHPParticleYield::GetDelayed(), G4ParticleHPParticleYield::GetMean(), G4ParticleHPParticleYield::GetPrompt(), and theFinalStateNeutrons.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

◆ SetA_Z()

void G4ParticleHPFinalState::SetA_Z ( G4double  anA,
G4double  aZ,
G4int  aM = 0 
)
inlineinherited

◆ SetAZMs()

void G4ParticleHPFinalState::SetAZMs ( G4double  anA,
G4double  aZ,
G4int  aM,
G4ParticleHPDataUsed  used 
)
inlineinherited

◆ SetNeutronRP()

void G4ParticleHPFSFissionFS::SetNeutronRP ( const G4ReactionProduct aNeutron)
inline

Definition at line 82 of file G4ParticleHPFSFissionFS.hh.

83 {
84 fCache.Get().theNeutronRP = &aNeutron;
86 }
void SetProjectileRP(const G4ReactionProduct &anIncidentParticleRP)

References fCache, G4ParticleHPAngular::SetProjectileRP(), and theNeutronAngularDis.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

◆ SetProjectile()

void G4ParticleHPFinalState::SetProjectile ( G4ParticleDefinition projectile)
inlineinherited

Definition at line 115 of file G4ParticleHPFinalState.hh.

116 {
117 theProjectile = projectile;
118 }

References G4ParticleHPFinalState::theProjectile.

Referenced by G4ParticleHPChannel::Register().

◆ SetTarget()

void G4ParticleHPFSFissionFS::SetTarget ( const G4ReactionProduct aTarget)
inline

Definition at line 88 of file G4ParticleHPFSFissionFS.hh.

89 {
90 fCache.Get().theTarget = &aTarget;
92 }
void SetTarget(const G4ReactionProduct &aTarget)

References fCache, G4ParticleHPAngular::SetTarget(), and theNeutronAngularDis.

Referenced by G4ParticleHPFissionFS::ApplyYourself().

Field Documentation

◆ fCache

G4Cache<toBeCached> G4ParticleHPFSFissionFS::fCache
private

Definition at line 113 of file G4ParticleHPFSFissionFS.hh.

Referenced by ApplyYourself(), GetPhotons(), SetNeutronRP(), and SetTarget().

◆ hasAnyData

G4bool G4ParticleHPFinalState::hasAnyData
protectedinherited

◆ hasFSData

G4bool G4ParticleHPFinalState::hasFSData
protectedinherited

◆ hasXsec

G4bool G4ParticleHPFinalState::hasXsec
protectedinherited

◆ secID

G4int G4ParticleHPFinalState::secID
protectedinherited

Definition at line 140 of file G4ParticleHPFinalState.hh.

Referenced by G4ParticleHPFinalState::adjust_final_state(), G4ParticleHPCaptureFS::ApplyYourself(), G4ParticleHPElasticFS::ApplyYourself(), G4ParticleHPFissionFS::ApplyYourself(), G4ParticleHPInelasticBaseFS::BaseApply(), G4ParticleHPInelasticCompFS::CompositeApply(), G4ParticleHP2AInelasticFS::G4ParticleHP2AInelasticFS(), G4ParticleHP2N2AInelasticFS::G4ParticleHP2N2AInelasticFS(), G4ParticleHP2NAInelasticFS::G4ParticleHP2NAInelasticFS(), G4ParticleHP2NDInelasticFS::G4ParticleHP2NDInelasticFS(), G4ParticleHP2NInelasticFS::G4ParticleHP2NInelasticFS(), G4ParticleHP2NPInelasticFS::G4ParticleHP2NPInelasticFS(), G4ParticleHP2PInelasticFS::G4ParticleHP2PInelasticFS(), G4ParticleHP3AInelasticFS::G4ParticleHP3AInelasticFS(), G4ParticleHP3NAInelasticFS::G4ParticleHP3NAInelasticFS(), G4ParticleHP3NInelasticFS::G4ParticleHP3NInelasticFS(), G4ParticleHP3NPInelasticFS::G4ParticleHP3NPInelasticFS(), G4ParticleHP4NInelasticFS::G4ParticleHP4NInelasticFS(), G4ParticleHPAInelasticFS::G4ParticleHPAInelasticFS(), G4ParticleHPCaptureFS::G4ParticleHPCaptureFS(), G4ParticleHPD2AInelasticFS::G4ParticleHPD2AInelasticFS(), G4ParticleHPDAInelasticFS::G4ParticleHPDAInelasticFS(), G4ParticleHPDInelasticFS::G4ParticleHPDInelasticFS(), G4ParticleHPElasticFS::G4ParticleHPElasticFS(), G4ParticleHPFinalState::G4ParticleHPFinalState(), G4ParticleHPFissionFS::G4ParticleHPFissionFS(), G4ParticleHPHe3InelasticFS::G4ParticleHPHe3InelasticFS(), G4ParticleHPN2AInelasticFS::G4ParticleHPN2AInelasticFS(), G4ParticleHPN2PInelasticFS::G4ParticleHPN2PInelasticFS(), G4ParticleHPN3AInelasticFS::G4ParticleHPN3AInelasticFS(), G4ParticleHPNAInelasticFS::G4ParticleHPNAInelasticFS(), G4ParticleHPND2AInelasticFS::G4ParticleHPND2AInelasticFS(), G4ParticleHPNDInelasticFS::G4ParticleHPNDInelasticFS(), G4ParticleHPNHe3InelasticFS::G4ParticleHPNHe3InelasticFS(), G4ParticleHPNInelasticFS::G4ParticleHPNInelasticFS(), G4ParticleHPNPAInelasticFS::G4ParticleHPNPAInelasticFS(), G4ParticleHPNPInelasticFS::G4ParticleHPNPInelasticFS(), G4ParticleHPNT2AInelasticFS::G4ParticleHPNT2AInelasticFS(), G4ParticleHPNTInelasticFS::G4ParticleHPNTInelasticFS(), G4ParticleHPNXInelasticFS::G4ParticleHPNXInelasticFS(), G4ParticleHPPAInelasticFS::G4ParticleHPPAInelasticFS(), G4ParticleHPPDInelasticFS::G4ParticleHPPDInelasticFS(), G4ParticleHPPInelasticFS::G4ParticleHPPInelasticFS(), G4ParticleHPPTInelasticFS::G4ParticleHPPTInelasticFS(), G4ParticleHPT2AInelasticFS::G4ParticleHPT2AInelasticFS(), G4ParticleHPTInelasticFS::G4ParticleHPTInelasticFS(), and G4ParticleHPInelasticCompFS::use_nresp71_model().

◆ theBaseA

G4double G4ParticleHPFinalState::theBaseA
protectedinherited

◆ theBaseM

G4int G4ParticleHPFinalState::theBaseM
protectedinherited

◆ theBaseZ

G4double G4ParticleHPFinalState::theBaseZ
protectedinherited

◆ theDelayedNeutronEnDis

G4ParticleHPEnergyDistribution G4ParticleHPFSFissionFS::theDelayedNeutronEnDis
private

Definition at line 107 of file G4ParticleHPFSFissionFS.hh.

Referenced by ApplyYourself(), and Init().

◆ theEnergyRelease

G4ParticleHPFissionERelease G4ParticleHPFSFissionFS::theEnergyRelease
private

Definition at line 111 of file G4ParticleHPFSFissionFS.hh.

Referenced by GetEnergyRelease(), and Init().

◆ theFinalStateNeutrons

G4ParticleHPParticleYield G4ParticleHPFSFissionFS::theFinalStateNeutrons
private

Definition at line 105 of file G4ParticleHPFSFissionFS.hh.

Referenced by ApplyYourself(), GetMass(), Init(), and SampleNeutronMult().

◆ theFinalStatePhotons

G4ParticleHPPhotonDist G4ParticleHPFSFissionFS::theFinalStatePhotons
private

Definition at line 110 of file G4ParticleHPFSFissionFS.hh.

Referenced by GetPhotons(), and Init().

◆ theNames

G4ParticleHPNames G4ParticleHPFSFissionFS::theNames
private

Definition at line 114 of file G4ParticleHPFSFissionFS.hh.

Referenced by Init().

◆ theNDLDataA

G4int G4ParticleHPFinalState::theNDLDataA
protectedinherited

◆ theNDLDataM

G4int G4ParticleHPFinalState::theNDLDataM
protectedinherited

◆ theNDLDataZ

G4int G4ParticleHPFinalState::theNDLDataZ
protectedinherited

◆ theNeutronAngularDis

G4ParticleHPAngular G4ParticleHPFSFissionFS::theNeutronAngularDis
private

Definition at line 108 of file G4ParticleHPFSFissionFS.hh.

Referenced by ApplyYourself(), Init(), SetNeutronRP(), and SetTarget().

◆ theProjectile

G4ParticleDefinition* G4ParticleHPFinalState::theProjectile
protectedinherited

◆ thePromptNeutronEnDis

G4ParticleHPEnergyDistribution G4ParticleHPFSFissionFS::thePromptNeutronEnDis
private

Definition at line 106 of file G4ParticleHPFSFissionFS.hh.

Referenced by ApplyYourself(), and Init().

◆ theResult

G4Cache< G4HadFinalState* > G4ParticleHPFinalState::theResult
protectedinherited

Definition at line 129 of file G4ParticleHPFinalState.hh.

Referenced by G4ParticleHPFinalState::adjust_final_state(), G4FissionLibrary::ApplyYourself(), G4ParticleHP2AInelasticFS::ApplyYourself(), G4ParticleHP2N2AInelasticFS::ApplyYourself(), G4ParticleHP2NAInelasticFS::ApplyYourself(), G4ParticleHP2NDInelasticFS::ApplyYourself(), G4ParticleHP2NInelasticFS::ApplyYourself(), G4ParticleHP2NPInelasticFS::ApplyYourself(), G4ParticleHP2PInelasticFS::ApplyYourself(), G4ParticleHP3AInelasticFS::ApplyYourself(), G4ParticleHP3NAInelasticFS::ApplyYourself(), G4ParticleHP3NInelasticFS::ApplyYourself(), G4ParticleHP3NPInelasticFS::ApplyYourself(), G4ParticleHP4NInelasticFS::ApplyYourself(), G4ParticleHPAInelasticFS::ApplyYourself(), G4ParticleHPCaptureFS::ApplyYourself(), G4ParticleHPD2AInelasticFS::ApplyYourself(), G4ParticleHPDAInelasticFS::ApplyYourself(), G4ParticleHPDInelasticFS::ApplyYourself(), G4ParticleHPElasticFS::ApplyYourself(), G4ParticleHPFissionFS::ApplyYourself(), G4ParticleHPHe3InelasticFS::ApplyYourself(), G4ParticleHPN2AInelasticFS::ApplyYourself(), G4ParticleHPN2PInelasticFS::ApplyYourself(), G4ParticleHPN3AInelasticFS::ApplyYourself(), G4ParticleHPNAInelasticFS::ApplyYourself(), G4ParticleHPND2AInelasticFS::ApplyYourself(), G4ParticleHPNDInelasticFS::ApplyYourself(), G4ParticleHPNHe3InelasticFS::ApplyYourself(), G4ParticleHPNInelasticFS::ApplyYourself(), G4ParticleHPNPAInelasticFS::ApplyYourself(), G4ParticleHPNPInelasticFS::ApplyYourself(), G4ParticleHPNT2AInelasticFS::ApplyYourself(), G4ParticleHPNTInelasticFS::ApplyYourself(), G4ParticleHPNXInelasticFS::ApplyYourself(), G4ParticleHPPAInelasticFS::ApplyYourself(), G4ParticleHPPDInelasticFS::ApplyYourself(), G4ParticleHPPInelasticFS::ApplyYourself(), G4ParticleHPPTInelasticFS::ApplyYourself(), G4ParticleHPT2AInelasticFS::ApplyYourself(), G4ParticleHPTInelasticFS::ApplyYourself(), G4ParticleHPInelasticBaseFS::BaseApply(), G4ParticleHPInelasticCompFS::CompositeApply(), G4ParticleHPFinalState::G4ParticleHPFinalState(), G4ParticleHPInelasticCompFS::use_nresp71_model(), and G4ParticleHPFinalState::~G4ParticleHPFinalState().


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