Geant4-11
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes
G4ParticleHPNInelasticFS Class Reference

#include <G4ParticleHPNInelasticFS.hh>

Inheritance diagram for G4ParticleHPNInelasticFS:
G4ParticleHPInelasticCompFS G4ParticleHPFinalState

Public Member Functions

G4HadFinalStateApplyYourself (const G4HadProjectile &theTrack)
 
void CompositeApply (const G4HadProjectile &theTrack, G4ParticleDefinition *aHadron)
 
 G4ParticleHPNInelasticFS ()
 
G4double GetA ()
 
G4int GetM ()
 
G4double GetN ()
 
virtual G4ParticleHPVectorGetXsec ()
 
virtual G4double GetXsec (G4double anEnergy)
 
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)
 
void InitDistributionInitialState (G4ReactionProduct &anIncidentPart, G4ReactionProduct &aTarget, G4int it)
 
void InitGammas (G4double AR, G4double ZR)
 
G4ParticleHPFinalStateNew ()
 
G4int SelectExitChannel (G4double eKinetic)
 
void SetA_Z (G4double anA, G4double aZ, G4int aM=0)
 
void SetAZMs (G4double anA, G4double aZ, G4int aM, G4ParticleHPDataUsed used)
 
void SetProjectile (G4ParticleDefinition *projectile)
 
 ~G4ParticleHPNInelasticFS ()
 

Protected Member Functions

void adjust_final_state (G4LorentzVector)
 

Protected Attributes

G4String gammaPath
 
G4bool hasAnyData
 
G4bool hasFSData
 
G4bool hasXsec
 
std::vector< G4intLR
 
std::vector< G4doubleQI
 
G4int secID
 
G4ParticleHPAngulartheAngularDistribution [51]
 
G4double theBaseA
 
G4int theBaseM
 
G4double theBaseZ
 
G4ParticleHPEnAngCorrelationtheEnergyAngData [51]
 
G4ParticleHPEnergyDistributiontheEnergyDistribution [51]
 
G4ParticleHPPhotonDisttheFinalStatePhotons [51]
 
G4ParticleHPDeExGammas theGammas
 
G4ParticleHPNames theNames
 
G4int theNDLDataA
 
G4int theNDLDataM
 
G4int theNDLDataZ
 
G4ParticleDefinitiontheProjectile
 
G4Cache< G4HadFinalState * > theResult
 
G4ParticleHPVectortheXsection [51]
 

Private Member Functions

void two_body_reaction (G4ReactionProduct *proj, G4ReactionProduct *targ, G4ReactionProduct *product, G4double resExcitationEnergy)
 
G4bool use_nresp71_model (const G4ParticleDefinition *aDefinition, const G4int it, const G4ReactionProduct &theTarget, G4ReactionProduct &boosted)
 

Private Attributes

G4NRESP71M03 nresp71_model
 

Detailed Description

Definition at line 41 of file G4ParticleHPNInelasticFS.hh.

Constructor & Destructor Documentation

◆ G4ParticleHPNInelasticFS()

G4ParticleHPNInelasticFS::G4ParticleHPNInelasticFS ( )

Definition at line 42 of file G4ParticleHPNInelasticFS.cc.

43{
44 secID = G4PhysicsModelCatalog::GetModelID( "model_G4ParticleHPNInelasticFS_F01" );
45}
static G4int GetModelID(const G4int modelIndex)

References G4PhysicsModelCatalog::GetModelID(), and G4ParticleHPFinalState::secID.

Referenced by New().

◆ ~G4ParticleHPNInelasticFS()

G4ParticleHPNInelasticFS::~G4ParticleHPNInelasticFS ( )
inline

Definition at line 46 of file G4ParticleHPNInelasticFS.hh.

46{}

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()

G4HadFinalState * G4ParticleHPNInelasticFS::ApplyYourself ( const G4HadProjectile theTrack)
virtual

Implements G4ParticleHPInelasticCompFS.

Definition at line 76 of file G4ParticleHPNInelasticFS.cc.

77{
78
79// do the final state
81
82// return the result
83 return theResult.Get();
84}
void CompositeApply(const G4HadProjectile &theTrack, G4ParticleDefinition *aHadron)

References G4ParticleHPInelasticCompFS::CompositeApply(), G4Cache< VALTYPE >::Get(), G4Neutron::Neutron(), and G4ParticleHPFinalState::theResult.

◆ CompositeApply()

void G4ParticleHPInelasticCompFS::CompositeApply ( const G4HadProjectile theTrack,
G4ParticleDefinition aHadron 
)
inherited

Definition at line 235 of file G4ParticleHPInelasticCompFS.cc.

237{
238
239// prepare neutron
240 if ( theResult.Get() == NULL ) theResult.Put( new G4HadFinalState );
241 theResult.Get()->Clear();
242 G4double eKinetic = theTrack.GetKineticEnergy();
243 const G4HadProjectile *hadProjectile = &theTrack;
244 G4ReactionProduct incidReactionProduct( const_cast<G4ParticleDefinition *>(hadProjectile->GetDefinition()) ); // incidReactionProduct
245 incidReactionProduct.SetMomentum( hadProjectile->Get4Momentum().vect() );
246 incidReactionProduct.SetKineticEnergy( eKinetic );
247
248// prepare target
249 G4int i;
250 for(i=0; i<50; i++)
251 { if(theXsection[i] != 0) { break; } }
252
253 G4double targetMass=0;
254 G4double eps = 0.0001;
255 targetMass = G4NucleiProperties::GetNuclearMass(static_cast<G4int>(theBaseA+eps), static_cast<G4int>(theBaseZ+eps));
256#ifdef G4PHPDEBUG
257 if( std::getenv("G4ParticleHPDebug")) G4cout <<this <<" G4ParticleHPInelasticCompFS::CompositeApply A " <<theBaseA <<" Z " <<theBaseZ <<" incident " <<hadProjectile->GetDefinition()->GetParticleName() <<G4endl;
258#endif
259// if(theEnergyAngData[i]!=0)
260// targetMass = theEnergyAngData[i]->GetTargetMass();
261// else if(theAngularDistribution[i]!=0)
262// targetMass = theAngularDistribution[i]->GetTargetMass();
263// else if(theFinalStatePhotons[50]!=0)
264// targetMass = theFinalStatePhotons[50]->GetTargetMass();
265 G4ReactionProduct theTarget;
266 G4Nucleus aNucleus;
267 //G4ThreeVector neuVelo = (1./hadProjectile->GetDefinition()->GetPDGMass())*incidReactionProduct.GetMomentum();
268 //theTarget = aNucleus.GetBiasedThermalNucleus( targetMass/hadProjectile->GetDefinition()->GetPDGMass() , neuVelo, theTrack.GetMaterial()->GetTemperature());
269 //G4Nucleus::GetBiasedThermalNucleus requests normalization of mass and velocity in neutron mass
270 G4ThreeVector neuVelo = ( 1./G4Neutron::Neutron()->GetPDGMass() )*incidReactionProduct.GetMomentum();
271 theTarget = aNucleus.GetBiasedThermalNucleus( targetMass/G4Neutron::Neutron()->GetPDGMass()
272 , neuVelo, theTrack.GetMaterial()->GetTemperature() );
273
275
276// prepare the residual mass
277 G4double residualMass=0;
278 G4double residualZ = theBaseZ + theProjectile->GetPDGCharge() - aDefinition->GetPDGCharge();
279 G4double residualA = theBaseA + theProjectile->GetBaryonNumber() - aDefinition->GetBaryonNumber();
280 residualMass = G4NucleiProperties::GetNuclearMass(static_cast<G4int>(residualA+eps), static_cast<G4int>(residualZ+eps));
281
282// prepare energy in target rest frame
283 G4ReactionProduct boosted;
284 boosted.Lorentz(incidReactionProduct, theTarget);
285 eKinetic = boosted.GetKineticEnergy();
286// G4double momentumInCMS = boosted.GetTotalMomentum();
287
288// select exit channel for composite FS class.
289 G4int it = SelectExitChannel( eKinetic );
290
291 //E. Mendoza (2018) -- to use JENDL/AN-2005
294 it=50;
295 }
296 }
297
298 // set target and neutron in the relevant exit channel
299 InitDistributionInitialState(incidReactionProduct, theTarget, it);
300
301 //---------------------------------------------------------------------//
302 //Hook for NRESP71MODEL
303 if ( G4ParticleHPManager::GetInstance()->GetUseNRESP71Model() && eKinetic<20*MeV) {
304 if ( (G4int)(theBaseZ+0.1) == 6 ) // If the reaction is with Carbon...
305 {
307 if ( use_nresp71_model( aDefinition , it , theTarget , boosted ) ) return;
308 }
309 }
310 }
311 //---------------------------------------------------------------------//
312
313 G4ReactionProductVector * thePhotons = 0;
314 G4ReactionProductVector * theParticles = 0;
315 G4ReactionProduct aHadron;
316 aHadron.SetDefinition(aDefinition); // what if only cross-sections exist ==> Na 23 11 @@@@
317 G4double availableEnergy = incidReactionProduct.GetKineticEnergy() + incidReactionProduct.GetMass() - aHadron.GetMass() +
318 (targetMass - residualMass);
319//080730c
320 if ( availableEnergy < 0 )
321 {
322 //G4cout << "080730c Adjust availavleEnergy " << G4endl;
323 availableEnergy = 0;
324 }
325 G4int nothingWasKnownOnHadron = 0;
326 G4int dummy;
327 G4double eGamm = 0;
328 G4int iLevel=it-1;
329
330// TK without photon has it = 0
331 if( 50 == it )
332 {
333
334// TK Excitation level is not determined
335 iLevel=-1;
336 aHadron.SetKineticEnergy(availableEnergy*residualMass/
337 (aHadron.GetMass()+residualMass));
338
339 //aHadron.SetMomentum(incidReactionProduct.GetMomentum()*(1./incidReactionProduct.GetTotalMomentum())*
340 // std::sqrt(aHadron.GetTotalEnergy()*aHadron.GetTotalEnergy()-
341 // aHadron.GetMass()*aHadron.GetMass()));
342
343 //TK add safty 100909
344 G4double p2 = ( aHadron.GetTotalEnergy()*aHadron.GetTotalEnergy() - aHadron.GetMass()*aHadron.GetMass() );
345 G4double p = 0.0;
346 if ( p2 > 0.0 ) p = std::sqrt( p );
347
348 aHadron.SetMomentum(incidReactionProduct.GetMomentum()*(1./incidReactionProduct.GetTotalMomentum())*p );
349
350 }
351 else
352 {
353 while ( iLevel!=-1 && theGammas.GetLevel(iLevel) == 0 ) { iLevel--; } // Loop checking, 11.05.2015, T. Koi
354 }
355
356
357 if ( theAngularDistribution[it] != 0 ) // MF4
358 {
359 if(theEnergyDistribution[it]!=0) // MF5
360 {
361 //************************************************************
362 /*
363 aHadron.SetKineticEnergy(theEnergyDistribution[it]->Sample(eKinetic, dummy));
364 G4double eSecN = aHadron.GetKineticEnergy();
365 */
366 //************************************************************
367 //EMendoza --> maximum allowable energy should be taken into account.
368 G4double dqi = 0.0;
369 if ( QI[it] < 0 || 849 < QI[it] ) dqi = QI[it]; //For backword compatibility QI introduced since G4NDL3.15
370 G4double MaxEne=eKinetic+dqi;
371 G4double eSecN=0.;
372
373 G4int icounter=0;
374 G4int icounter_max=1024;
375 do {
376 icounter++;
377 if ( icounter > icounter_max ) {
378 G4cout << "Loop-counter exceeded the threshold value at " << __LINE__ << "th line of " << __FILE__ << "." << G4endl;
379 break;
380 }
381 eSecN=theEnergyDistribution[it]->Sample(eKinetic, dummy);
382 }while(eSecN>MaxEne); // Loop checking, 11.05.2015, T. Koi
383 aHadron.SetKineticEnergy(eSecN);
384 //************************************************************
385 eGamm = eKinetic-eSecN;
386 for(iLevel=theGammas.GetNumberOfLevels()-1; iLevel>=0; iLevel--)
387 {
388 if(theGammas.GetLevelEnergy(iLevel)<eGamm) break;
389 }
390 G4double random = 2*G4UniformRand();
391 iLevel+=G4int(random);
392 if(iLevel>theGammas.GetNumberOfLevels()-1)iLevel = theGammas.GetNumberOfLevels()-1;
393 }
394 else
395 {
396 G4double eExcitation = 0;
397 if(iLevel>=0) eExcitation = theGammas.GetLevel(iLevel)->GetLevelEnergy();
398 while (eKinetic-eExcitation < 0 && iLevel>0) // Loop checking, 11.05.2015, T. Koi
399 {
400 iLevel--;
401 eExcitation = theGammas.GetLevel(iLevel)->GetLevelEnergy();
402 }
403 //110610TK BEGIN
404 //Use QI value for calculating excitation energy of residual.
405 G4bool useQI=false;
406 G4double dqi = QI[it];
407 if (dqi < 0 || 849 < dqi) useQI = true; // Former libraries do not have values in this range
408
409 if (useQI) {
410 // QI introudced since G4NDL3.15
411 // G4double QM=(incidReactionProduct.GetMass()+targetMass)-(aHadron.GetMass()+residualMass);
412 // eExcitation = QM-QI[it];
413 // eExcitation = QI[0] - QI[it]; // Bug fix #1838
414 // if(eExcitation < 20*CLHEP::keV) eExcitation = 0;
415
416 eExcitation = std::max(0.,QI[0] - QI[it]); // Bug fix 2333
417
418 // Re-evluate iLevel based on this eExcitation
419 iLevel = 0;
420 G4bool find = false;
421 G4int imaxEx = 0;
422 G4double level_tolerance = 1.0*CLHEP::keV;
423
424 while( theGammas.GetLevel(iLevel+1) != 0 ) { // Loop checking, 11.05.2015, T. Koi
425 G4double maxEx = 0.0;
426 if (maxEx < theGammas.GetLevel(iLevel)->GetLevelEnergy() ) {
427 maxEx = theGammas.GetLevel(iLevel)->GetLevelEnergy();
428 imaxEx = iLevel;
429 }
430
431 // Fix bug 1789 DHW - first if-branch added because gamma data come from ENSDF
432 // and do not necessarily match the excitations used in ENDF-B.VII
433 // Compromise solution: use 1 keV tolerance suggested by T. Koi
434 if (std::abs(eExcitation - theGammas.GetLevel(iLevel)->GetLevelEnergy() ) < level_tolerance) {
435 find = true;
436 break;
437
438 } else if (eExcitation < theGammas.GetLevel(iLevel)->GetLevelEnergy() ) {
439 find = true;
440 iLevel--;
441 // very small eExcitation, iLevel becomes -1, this is protected below
442 if (theTrack.GetDefinition() == aDefinition) { // this line added as part of fix #1838
443 if (iLevel == -1) iLevel = 0;
444 }
445 break;
446 }
447 iLevel++;
448 }
449
450 // If proper level cannot be found, use the maximum level
451 if (!find) iLevel = imaxEx;
452 }
453
454 if(std::getenv("G4ParticleHPDebug") && eKinetic-eExcitation < 0)
455 {
456 throw G4HadronicException(__FILE__, __LINE__, "SEVERE: InelasticCompFS: Consistency of data not good enough, please file report");
457 }
458 if(eKinetic-eExcitation < 0) eExcitation = 0;
459 if(iLevel!= -1) aHadron.SetKineticEnergy(eKinetic - eExcitation);
460
461 }
463
464 if( theFinalStatePhotons[it] == 0 )
465 {
466 //G4cout << "110610 USE Gamma Level" << G4endl;
467// TK comment Most n,n* eneter to this
468 thePhotons = theGammas.GetDecayGammas(iLevel);
469 eGamm -= theGammas.GetLevelEnergy(iLevel);
470 if(eGamm>0) // @ ok for now, but really needs an efficient way of correllated sampling @
471 {
472 G4ReactionProduct * theRestEnergy = new G4ReactionProduct;
473 theRestEnergy->SetDefinition(G4Gamma::Gamma());
474 theRestEnergy->SetKineticEnergy(eGamm);
475 G4double costh = 2.*G4UniformRand()-1.;
477 theRestEnergy->SetMomentum(eGamm*std::sin(std::acos(costh))*std::cos(phi),
478 eGamm*std::sin(std::acos(costh))*std::sin(phi),
479 eGamm*costh);
480 if(thePhotons == 0) { thePhotons = new G4ReactionProductVector; }
481 thePhotons->push_back(theRestEnergy);
482 }
483 }
484 }
485 else if(theEnergyAngData[it] != 0) // MF6
486 {
487
488 theParticles = theEnergyAngData[it]->Sample(eKinetic);
489
490 //141017 Fix BEGIN
491 //Adjust A and Z in the case of miss much between selected data and target nucleus
492 if ( theParticles != NULL ) {
493 G4int sumA = 0;
494 G4int sumZ = 0;
495 G4int maxA = 0;
496 G4int jAtMaxA = 0;
497 for ( G4int j = 0 ; j != (G4int)theParticles->size() ; j++ ) {
498 if ( theParticles->at(j)->GetDefinition()->GetBaryonNumber() > maxA ) {
499 maxA = theParticles->at(j)->GetDefinition()->GetBaryonNumber();
500 jAtMaxA = j;
501 }
502 sumA += theParticles->at(j)->GetDefinition()->GetBaryonNumber();
503 sumZ += G4int( theParticles->at(j)->GetDefinition()->GetPDGCharge() + eps );
504 }
505 G4int dA = (G4int)theBaseA + hadProjectile->GetDefinition()->GetBaryonNumber() - sumA;
506 G4int dZ = (G4int)theBaseZ + G4int( hadProjectile->GetDefinition()->GetPDGCharge() + eps ) - sumZ;
507 if ( dA < 0 || dZ < 0 ) {
508 G4int newA = theParticles->at(jAtMaxA)->GetDefinition()->GetBaryonNumber() + dA ;
509 G4int newZ = G4int( theParticles->at(jAtMaxA)->GetDefinition()->GetPDGCharge() + eps ) + dZ;
510 G4ParticleDefinition* pd = G4IonTable::GetIonTable()->GetIon ( newZ , newA );
511 theParticles->at( jAtMaxA )->SetDefinition( pd );
512 }
513 }
514 //141017 Fix END
515
516 }
517 else
518 {
519 // @@@ what to do, if we have photon data, but no info on the hadron itself
520 nothingWasKnownOnHadron = 1;
521 }
522
523 //G4cout << "theFinalStatePhotons it " << it << G4endl;
524 //G4cout << "theFinalStatePhotons[it] " << theFinalStatePhotons[it] << G4endl;
525 //G4cout << "theFinalStatePhotons it " << it << G4endl;
526 //G4cout << "theFinalStatePhotons[it] " << theFinalStatePhotons[it] << G4endl;
527 //G4cout << "thePhotons " << thePhotons << G4endl;
528
529 if ( theFinalStatePhotons[it] != 0 )
530 {
531 // the photon distributions are in the Nucleus rest frame.
532 // TK residual rest frame
533 G4ReactionProduct boosted_tmp;
534 boosted_tmp.Lorentz(incidReactionProduct, theTarget);
535 G4double anEnergy = boosted_tmp.GetKineticEnergy();
536 thePhotons = theFinalStatePhotons[it]->GetPhotons(anEnergy);
537 G4double aBaseEnergy = theFinalStatePhotons[it]->GetLevelEnergy();
538 G4double testEnergy = 0;
539 if(thePhotons!=0 && thePhotons->size()!=0)
540 { aBaseEnergy-=thePhotons->operator[](0)->GetTotalEnergy(); }
541 if(theFinalStatePhotons[it]->NeedsCascade())
542 {
543 while(aBaseEnergy>0.01*CLHEP::keV) // Loop checking, 11.05.2015, T. Koi
544 {
545 // cascade down the levels
546 G4bool foundMatchingLevel = false;
547 G4int closest = 2;
548 G4double deltaEold = -1;
549 for(G4int j=1; j<it; j++)
550 {
551 if(theFinalStatePhotons[j]!=0)
552 {
553 testEnergy = theFinalStatePhotons[j]->GetLevelEnergy();
554 }
555 else
556 {
557 testEnergy = 0;
558 }
559 G4double deltaE = std::abs(testEnergy-aBaseEnergy);
560 if(deltaE<0.1*CLHEP::keV)
561 {
562 G4ReactionProductVector * theNext =
563 theFinalStatePhotons[j]->GetPhotons(anEnergy);
564 if ( thePhotons != NULL ) thePhotons->push_back(theNext->operator[](0));
565 aBaseEnergy = testEnergy-theNext->operator[](0)->GetTotalEnergy();
566 delete theNext;
567 foundMatchingLevel = true;
568 break; // ===>
569 }
570 if(theFinalStatePhotons[j]!=0 && ( deltaE<deltaEold||deltaEold<0.) )
571 {
572 closest = j;
573 deltaEold = deltaE;
574 }
575 } // <=== the break goes here.
576 if(!foundMatchingLevel)
577 {
578 G4ReactionProductVector * theNext =
579 theFinalStatePhotons[closest]->GetPhotons(anEnergy);
580 if ( thePhotons != NULL ) thePhotons->push_back(theNext->operator[](0));
581 aBaseEnergy = aBaseEnergy-theNext->operator[](0)->GetTotalEnergy();
582 delete theNext;
583 }
584 }
585 }
586 }
587 unsigned int i0;
588 if(thePhotons!=0)
589 {
590 for(i0=0; i0<thePhotons->size(); i0++)
591 {
592 // back to lab
593 thePhotons->operator[](i0)->Lorentz(*(thePhotons->operator[](i0)), -1.*theTarget);
594 }
595 }
596 //G4cout << "nothingWasKnownOnHadron " << nothingWasKnownOnHadron << G4endl;
597 if (nothingWasKnownOnHadron)
598 {
599// In this case, hadron should be isotropic in CM
600// Next 12 lines are Emilio's replacement
601 // G4double QM=(incidReactionProduct.GetMass()+targetMass)-(aHadron.GetMass()+residualMass);
602 // G4double eExcitation = QM-QI[it];
603 // G4double eExcitation = QI[0] - QI[it]; // Fix of bug #1838
604 // if(eExcitation<20*CLHEP::keV){eExcitation=0;}
605
606 G4double eExcitation = std::max(0.,QI[0] - QI[it]); // Fix of bug #2333
607
608 two_body_reaction(&incidReactionProduct,&theTarget,&aHadron,eExcitation);
609 if(thePhotons==0 && eExcitation>0){
610 for(iLevel=theGammas.GetNumberOfLevels()-1; iLevel>=0; iLevel--)
611 {
612 if(theGammas.GetLevelEnergy(iLevel)<eExcitation+5*keV) break; // 5 keV tolerance
613 }
614 thePhotons = theGammas.GetDecayGammas(iLevel);
615 }
616 }
617// Emilio's replacement done
618/*
619// This code replaced by Emilio (previous 12 lines)
620// mu and p should be correlated
621//
622 //isotropic distribution in CM
623 G4double mu = 1.0 - 2.*G4UniformRand();
624
625 // Need momenta in target rest frame
626 G4LorentzVector target_in_LAB ( theTarget.GetMomentum() , theTarget.GetTotalEnergy() );
627 G4ThreeVector boostToTargetRest = -target_in_LAB.boostVector();
628 G4LorentzVector proj_in_LAB = hadProjectile->Get4Momentum();
629
630 G4DynamicParticle* proj = new G4DynamicParticle(theProjectile, proj_in_LAB.boost(boostToTargetRest) );
631// G4DynamicParticle* targ =
632// new G4DynamicParticle(G4IonTable::GetIonTable()->GetIon((G4int)theBaseZ, (G4int)theBaseA, totalPhotonEnergy), G4ThreeVector(0) );
633// Fix bug 2166 (A. Zontikov): replace above two lines with next three lines
634 G4double excitationEnergy = theFinalStatePhotons[it] ? theFinalStatePhotons[it]->GetLevelEnergy() : 0.0;
635 G4DynamicParticle* targ =
636 new G4DynamicParticle(G4IonTable::GetIonTable()->GetIon((G4int)theBaseZ, (G4int)theBaseA, excitationEnergy), G4ThreeVector(0) );
637 G4DynamicParticle* hadron =
638 new G4DynamicParticle(aHadron.GetDefinition(), G4ThreeVector(0) ); // Will fill in the momentum
639
640 two_body_reaction ( proj , targ , hadron , mu );
641
642 G4LorentzVector hadron_in_trag_rest = hadron->Get4Momentum();
643 G4LorentzVector hadron_in_LAB = hadron_in_trag_rest.boost ( -boostToTargetRest );
644 aHadron.SetMomentum( hadron_in_LAB.v() );
645 aHadron.SetKineticEnergy ( hadron_in_LAB.e() - hadron_in_LAB.m() );
646
647 delete proj;
648 delete targ;
649 delete hadron;
650
651 }
652*/
653
654// fill the result
655// Beware - the recoil is not necessarily in the particles...
656// Can be calculated from momentum conservation?
657// The idea is that the particles ar emitted forst, and the gammas only once the
658// recoil is on the residual; assumption is that gammas do not contribute to
659// the recoil.
660// This needs more design @@@
661
662 G4int nSecondaries = 2; // the hadron and the recoil
663 G4bool needsSeparateRecoil = false;
664 G4int totalBaryonNumber = 0;
665 G4int totalCharge = 0;
666 G4ThreeVector totalMomentum(0);
667 if(theParticles != 0)
668 {
669 nSecondaries = theParticles->size();
670 const G4ParticleDefinition * aDef;
671 unsigned int ii0;
672 for(ii0=0; ii0<theParticles->size(); ii0++)
673 {
674 aDef = theParticles->operator[](ii0)->GetDefinition();
675 totalBaryonNumber+=aDef->GetBaryonNumber();
676 totalCharge+=G4int(aDef->GetPDGCharge()+eps);
677 totalMomentum += theParticles->operator[](ii0)->GetMomentum();
678 }
679 if(totalBaryonNumber!=G4int(theBaseA+eps+hadProjectile->GetDefinition()->GetBaryonNumber()))
680 {
681 needsSeparateRecoil = true;
682 nSecondaries++;
683 residualA = G4int(theBaseA+eps+hadProjectile->GetDefinition()->GetBaryonNumber()
684 -totalBaryonNumber);
685 residualZ = G4int(theBaseZ+eps+hadProjectile->GetDefinition()->GetPDGCharge()
686 -totalCharge);
687 }
688 }
689
690 G4int nPhotons = 0;
691 if(thePhotons!=0) { nPhotons = thePhotons->size(); }
692 nSecondaries += nPhotons;
693
694 G4DynamicParticle * theSec;
695
696 if( theParticles==0 )
697 {
698 theSec = new G4DynamicParticle;
699 theSec->SetDefinition(aHadron.GetDefinition());
700 theSec->SetMomentum(aHadron.GetMomentum());
701 theResult.Get()->AddSecondary(theSec, secID);
702#ifdef G4PHPDEBUG
703 if( std::getenv("G4ParticleHPDebug")) G4cout << this << " G4ParticleHPInelasticCompFS::BaseApply add secondary1 " << theSec->GetParticleDefinition()->GetParticleName() << " E= " << theSec->GetKineticEnergy() << " NSECO " << theResult.Get()->GetNumberOfSecondaries() << G4endl;
704#endif
705
706 aHadron.Lorentz(aHadron, theTarget);
707 G4ReactionProduct theResidual;
709 ->GetIon(static_cast<G4int>(residualZ), static_cast<G4int>(residualA), 0));
710 theResidual.SetKineticEnergy(aHadron.GetKineticEnergy()*aHadron.GetMass()/theResidual.GetMass());
711
712 //080612TK contribution from Benoit Pirard and Laurent Desorgher (Univ. Bern) #6
713 //theResidual.SetMomentum(-1.*aHadron.GetMomentum());
714 G4ThreeVector incidentNeutronMomentum = incidReactionProduct.GetMomentum();
715 theResidual.SetMomentum(incidentNeutronMomentum - aHadron.GetMomentum());
716
717 theResidual.Lorentz(theResidual, -1.*theTarget);
718 G4ThreeVector totalPhotonMomentum(0,0,0);
719 if(thePhotons!=0)
720 {
721 for(i=0; i<nPhotons; i++)
722 {
723 totalPhotonMomentum += thePhotons->operator[](i)->GetMomentum();
724 }
725 }
726 theSec = new G4DynamicParticle;
727 theSec->SetDefinition(theResidual.GetDefinition());
728 theSec->SetMomentum(theResidual.GetMomentum()-totalPhotonMomentum);
729 theResult.Get()->AddSecondary(theSec, secID);
730#ifdef G4PHPDEBUG
731 if( std::getenv("G4ParticleHPDebug")) G4cout << this << " G4ParticleHPInelasticCompFS::BaseApply add secondary2 " << theSec->GetParticleDefinition()->GetParticleName() << " E= " << theSec->GetKineticEnergy() << " NSECO " << theResult.Get()->GetNumberOfSecondaries() << G4endl;
732#endif
733 }
734 else
735 {
736 for(i0=0; i0<theParticles->size(); i0++)
737 {
738 theSec = new G4DynamicParticle;
739 theSec->SetDefinition(theParticles->operator[](i0)->GetDefinition());
740 theSec->SetMomentum(theParticles->operator[](i0)->GetMomentum());
741 theResult.Get()->AddSecondary(theSec, secID);
742#ifdef G4PHPDEBUG
743 if( std::getenv("G4ParticleHPDebug")) G4cout << this << " G4ParticleHPInelasticCompFS::BaseApply add secondary3 " << theSec->GetParticleDefinition()->GetParticleName() << " E= " << theSec->GetKineticEnergy() << " NSECO " << theResult.Get()->GetNumberOfSecondaries() << G4endl;
744#endif
745 delete theParticles->operator[](i0);
746 }
747 delete theParticles;
748 if(needsSeparateRecoil && residualZ!=0)
749 {
750 G4ReactionProduct theResidual;
752 ->GetIon(static_cast<G4int>(residualZ), static_cast<G4int>(residualA), 0));
753 G4double resiualKineticEnergy = theResidual.GetMass()*theResidual.GetMass();
754 resiualKineticEnergy += totalMomentum*totalMomentum;
755 resiualKineticEnergy = std::sqrt(resiualKineticEnergy) - theResidual.GetMass();
756// cout << "Kinetic energy of the residual = "<<resiualKineticEnergy<<endl;
757 theResidual.SetKineticEnergy(resiualKineticEnergy);
758
759 //080612TK contribution from Benoit Pirard and Laurent Desorgher (Univ. Bern) #4
760 //theResidual.SetMomentum(-1.*totalMomentum);
761 //G4ThreeVector incidentNeutronMomentum = incidReactionProduct.GetMomentum();
762 //theResidual.SetMomentum(incidentNeutronMomentum - aHadron.GetMomentum());
763//080717 TK Comment still do NOT include photon's mometum which produce by thePhotons
764 theResidual.SetMomentum( incidReactionProduct.GetMomentum() + theTarget.GetMomentum() - totalMomentum );
765
766 theSec = new G4DynamicParticle;
767 theSec->SetDefinition(theResidual.GetDefinition());
768 theSec->SetMomentum(theResidual.GetMomentum());
769 theResult.Get()->AddSecondary(theSec, secID);
770#ifdef G4PHPDEBUG
771 if( std::getenv("G4ParticleHPDebug")) G4cout << this << " G4ParticleHPInelasticCompFS::BaseApply add secondary4 " << theSec->GetParticleDefinition()->GetParticleName() << " E= " << theSec->GetKineticEnergy() << " NSECO " << theResult.Get()->GetNumberOfSecondaries() << G4endl;
772#endif
773
774 }
775 }
776 if(thePhotons!=0)
777 {
778 for(i=0; i<nPhotons; i++)
779 {
780 theSec = new G4DynamicParticle;
781 //Bug reported Chao Zhang (Chao.Zhang@usd.edu), Dongming Mei(Dongming.Mei@usd.edu) Feb. 25, 2009
782 //theSec->SetDefinition(G4Gamma::Gamma());
783 theSec->SetDefinition( thePhotons->operator[](i)->GetDefinition() );
784 //But never cause real effect at least with G4NDL3.13 TK
785 theSec->SetMomentum(thePhotons->operator[](i)->GetMomentum());
786 theResult.Get()->AddSecondary(theSec, secID);
787#ifdef G4PHPDEBUG
788 if( std::getenv("G4ParticleHPDebug")) G4cout << this << " G4ParticleHPInelasticCompFS::BaseApply add secondary5 " << theSec->GetParticleDefinition()->GetParticleName() << " E= " << theSec->GetKineticEnergy() << " NSECO " << theResult.Get()->GetNumberOfSecondaries() << G4endl;
789#endif
790
791 delete thePhotons->operator[](i);
792 }
793// some garbage collection
794 delete thePhotons;
795 }
796
797//080721
799 G4LorentzVector targ_4p_lab ( theTarget.GetMomentum() , std::sqrt( targ_pd->GetPDGMass()*targ_pd->GetPDGMass() + theTarget.GetMomentum().mag2() ) );
800 G4LorentzVector proj_4p_lab = theTrack.Get4Momentum();
801 G4LorentzVector init_4p_lab = proj_4p_lab + targ_4p_lab;
802 adjust_final_state ( init_4p_lab );
803
804// clean up the primary neutron
806}
static const G4double eps
static const G4int maxA
@ stopAndKill
std::vector< G4ReactionProduct * > G4ReactionProductVector
double mag2() const
Hep3Vector vect() const
void Put(const value_type &val) const
Definition: G4Cache.hh:321
const G4ParticleDefinition * GetParticleDefinition() const
void SetStatusChange(G4HadFinalStateStatus aS)
const G4Material * GetMaterial() const
const G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
G4double GetTemperature() const
Definition: G4Material.hh:178
static G4Neutron * Definition()
Definition: G4Neutron.cc:53
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4ReactionProduct GetBiasedThermalNucleus(G4double aMass, G4ThreeVector aVelocity, G4double temp=-1) const
Definition: G4Nucleus.cc:118
G4double GetPDGCharge() const
void SampleAndUpdate(G4ReactionProduct &anIncidentParticle)
G4ParticleHPLevel * GetLevel(G4int i)
G4ReactionProductVector * GetDecayGammas(G4int aLevel)
G4double GetLevelEnergy(G4int aLevel)
G4ReactionProductVector * Sample(G4double anEnergy)
G4double Sample(G4double anEnergy, G4int &it)
void adjust_final_state(G4LorentzVector)
void InitDistributionInitialState(G4ReactionProduct &anIncidentPart, G4ReactionProduct &aTarget, G4int it)
G4int SelectExitChannel(G4double eKinetic)
G4ParticleHPAngular * theAngularDistribution[51]
G4ParticleHPEnergyDistribution * theEnergyDistribution[51]
G4ParticleHPEnAngCorrelation * theEnergyAngData[51]
G4ParticleHPPhotonDist * theFinalStatePhotons[51]
void two_body_reaction(G4ReactionProduct *proj, G4ReactionProduct *targ, G4ReactionProduct *product, G4double resExcitationEnergy)
G4bool use_nresp71_model(const G4ParticleDefinition *aDefinition, const G4int it, const G4ReactionProduct &theTarget, G4ReactionProduct &boosted)
G4ReactionProductVector * GetPhotons(G4double anEnergy)
void SetMomentum(const G4double x, const G4double y, const G4double z)
G4double GetKineticEnergy() const
const G4ParticleDefinition * GetDefinition() const
G4double GetTotalEnergy() const
G4ThreeVector GetMomentum() const
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
void SetKineticEnergy(const G4double en)
G4double GetMass() const
static constexpr double keV
static constexpr double twopi
Definition: SystemOfUnits.h:56

References G4HadFinalState::AddSecondary(), G4ParticleHPFinalState::adjust_final_state(), G4HadFinalState::Clear(), G4Neutron::Definition(), eps, G4cout, G4endl, G4UniformRand, G4Gamma::Gamma(), G4Cache< VALTYPE >::Get(), G4HadProjectile::Get4Momentum(), G4ParticleDefinition::GetBaryonNumber(), G4Nucleus::GetBiasedThermalNucleus(), G4ParticleHPDeExGammas::GetDecayGammas(), G4HadProjectile::GetDefinition(), G4ReactionProduct::GetDefinition(), G4ParticleHPManager::GetInstance(), G4IonTable::GetIon(), G4IonTable::GetIonTable(), G4DynamicParticle::GetKineticEnergy(), G4HadProjectile::GetKineticEnergy(), G4ReactionProduct::GetKineticEnergy(), G4ParticleHPDeExGammas::GetLevel(), G4ParticleHPLevel::GetLevelEnergy(), G4ParticleHPPhotonDist::GetLevelEnergy(), G4ParticleHPDeExGammas::GetLevelEnergy(), G4ReactionProduct::GetMass(), G4HadProjectile::GetMaterial(), G4ReactionProduct::GetMomentum(), G4NucleiProperties::GetNuclearMass(), G4ParticleHPDeExGammas::GetNumberOfLevels(), G4HadFinalState::GetNumberOfSecondaries(), G4DynamicParticle::GetParticleDefinition(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), G4ParticleHPPhotonDist::GetPhotons(), G4Material::GetTemperature(), G4ReactionProduct::GetTotalEnergy(), G4ReactionProduct::GetTotalMomentum(), G4ParticleHPInelasticCompFS::InitDistributionInitialState(), CLHEP::keV, keV, G4ReactionProduct::Lorentz(), CLHEP::Hep3Vector::mag2(), G4INCL::Math::max(), maxA, MeV, G4Neutron::Neutron(), G4Cache< VALTYPE >::Put(), G4ParticleHPInelasticCompFS::QI, G4ParticleHPEnAngCorrelation::Sample(), G4ParticleHPEnergyDistribution::Sample(), G4ParticleHPAngular::SampleAndUpdate(), G4ParticleHPFinalState::secID, G4ParticleHPInelasticCompFS::SelectExitChannel(), G4DynamicParticle::SetDefinition(), G4ReactionProduct::SetDefinition(), G4ReactionProduct::SetKineticEnergy(), G4ReactionProduct::SetMomentum(), G4DynamicParticle::SetMomentum(), G4HadFinalState::SetStatusChange(), stopAndKill, G4ParticleHPInelasticCompFS::theAngularDistribution, G4ParticleHPFinalState::theBaseA, G4ParticleHPFinalState::theBaseZ, G4ParticleHPInelasticCompFS::theEnergyAngData, G4ParticleHPInelasticCompFS::theEnergyDistribution, G4ParticleHPInelasticCompFS::theFinalStatePhotons, G4ParticleHPInelasticCompFS::theGammas, G4ParticleHPFinalState::theProjectile, G4ParticleHPFinalState::theResult, G4ParticleHPInelasticCompFS::theXsection, G4ParticleHPInelasticCompFS::two_body_reaction(), CLHEP::twopi, G4ParticleHPInelasticCompFS::use_nresp71_model(), and CLHEP::HepLorentzVector::vect().

Referenced by G4ParticleHPAInelasticFS::ApplyYourself(), G4ParticleHPDInelasticFS::ApplyYourself(), G4ParticleHPHe3InelasticFS::ApplyYourself(), ApplyYourself(), G4ParticleHPPInelasticFS::ApplyYourself(), and G4ParticleHPTInelasticFS::ApplyYourself().

◆ GetA()

G4double G4ParticleHPFinalState::GetA ( void  )
inlineinherited

Definition at line 104 of file G4ParticleHPFinalState.hh.

104{ return theBaseA; }

References G4ParticleHPFinalState::theBaseA.

◆ GetM()

G4int G4ParticleHPFinalState::GetM ( )
inlineinherited

◆ GetN()

G4double G4ParticleHPFinalState::GetN ( )
inlineinherited

◆ GetXsec() [1/2]

virtual G4ParticleHPVector * G4ParticleHPInelasticCompFS::GetXsec ( )
inlinevirtualinherited

◆ GetXsec() [2/2]

virtual G4double G4ParticleHPInelasticCompFS::GetXsec ( G4double  anEnergy)
inlinevirtualinherited

Reimplemented from G4ParticleHPFinalState.

Definition at line 90 of file G4ParticleHPInelasticCompFS.hh.

91 {
92 return std::max(0., theXsection[50]->GetY(anEnergy));
93 }

References G4INCL::Math::max(), and G4ParticleHPInelasticCompFS::theXsection.

◆ GetZ()

G4double G4ParticleHPFinalState::GetZ ( void  )
inlineinherited

◆ HasAnyData()

G4bool G4ParticleHPFinalState::HasAnyData ( )
inlineinherited

◆ HasFSData()

G4bool G4ParticleHPFinalState::HasFSData ( )
inlineinherited

◆ HasXsec()

G4bool G4ParticleHPFinalState::HasXsec ( )
inlineinherited

◆ Init() [1/2]

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

Reimplemented from G4ParticleHPInelasticCompFS.

Definition at line 47 of file G4ParticleHPNInelasticFS.cc.

48{
49 // G4cout << "Getting initialized for: "<<A<<" "<< Z<<G4endl;
50 G4ParticleHPInelasticCompFS::Init(A, Z, M, dirName, aFSType, projectile);
51 G4double ResidualA = 0;
52 G4double ResidualZ = 0;
53 if( projectile == G4Neutron::Neutron() ) {
54 ResidualA = A;
55 ResidualZ = Z;
56 } else if( projectile == G4Proton::Proton() ) {
57 ResidualA = A;
58 ResidualZ = Z+1;
59 } else if( projectile == G4Deuteron::Deuteron() ) {
60 ResidualA = A+1;
61 ResidualZ = Z+1;
62 } else if( projectile == G4Triton::Triton() ) {
63 ResidualA = A+2;
64 ResidualZ = Z+1;
65 } else if( projectile == G4He3::He3() ) {
66 ResidualA = A+2;
67 ResidualZ = Z+2;
68 } else if( projectile == G4Alpha::Alpha() ) {
69 ResidualA = A+3;
70 ResidualZ = Z+2;
71 }
72
73 G4ParticleHPInelasticCompFS::InitGammas(ResidualA, ResidualZ);
74}
#define M(row, col)
const G4int Z[17]
const G4double A[17]
void InitGammas(G4double AR, G4double ZR)
void Init(G4double A, G4double Z, G4int M, G4String &dirName, G4String &aSFType, G4ParticleDefinition *)

References A, G4Alpha::Alpha(), G4Deuteron::Deuteron(), G4He3::He3(), G4ParticleHPInelasticCompFS::Init(), G4ParticleHPInelasticCompFS::InitGammas(), M, G4Neutron::Neutron(), G4Proton::Proton(), G4Triton::Triton(), and Z.

◆ 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().

◆ InitDistributionInitialState()

void G4ParticleHPInelasticCompFS::InitDistributionInitialState ( G4ReactionProduct anIncidentPart,
G4ReactionProduct aTarget,
G4int  it 
)
inlineinherited

◆ InitGammas()

void G4ParticleHPInelasticCompFS::InitGammas ( G4double  AR,
G4double  ZR 
)
inherited

Definition at line 64 of file G4ParticleHPInelasticCompFS.cc.

65{
66 // char the[100] = {""};
67 // std::ostrstream ost(the, 100, std::ios::out);
68 // ost <<gammaPath<<"z"<<ZR<<".a"<<AR;
69 // G4String * aName = new G4String(the);
70 // std::ifstream from(*aName, std::ios::in);
71
72 std::ostringstream ost;
73 ost <<gammaPath<<"z"<<ZR<<".a"<<AR;
74 G4String aName = ost.str();
75 std::ifstream from(aName, std::ios::in);
76
77 if(!from) return; // no data found for this isotope
78 // std::ifstream theGammaData(*aName, std::ios::in);
79 std::ifstream theGammaData(aName, std::ios::in);
80
81 theGammas.Init(theGammaData);
82 // delete aName;
83
84}
void Init(std::istream &aDataFile)

References G4ParticleHPInelasticCompFS::gammaPath, G4ParticleHPDeExGammas::Init(), and G4ParticleHPInelasticCompFS::theGammas.

Referenced by G4ParticleHPAInelasticFS::Init(), G4ParticleHPDInelasticFS::Init(), G4ParticleHPHe3InelasticFS::Init(), Init(), G4ParticleHPPInelasticFS::Init(), and G4ParticleHPTInelasticFS::Init().

◆ New()

G4ParticleHPFinalState * G4ParticleHPNInelasticFS::New ( )
inlinevirtual

◆ SelectExitChannel()

G4int G4ParticleHPInelasticCompFS::SelectExitChannel ( G4double  eKinetic)
inherited

Definition at line 200 of file G4ParticleHPInelasticCompFS.cc.

201{
202
203// it = 0 has without Photon
204 G4double running[50];
205 running[0] = 0;
206 unsigned int i;
207 for(i=0; i<50; i++)
208 {
209 if(i!=0) running[i]=running[i-1];
210 if(theXsection[i] != 0)
211 {
212 running[i] += std::max(0., theXsection[i]->GetXsec(eKinetic));
213 }
214 }
215 G4double random = G4UniformRand();
216 G4double sum = running[49];
217 G4int it = 50;
218 if(0!=sum)
219 {
220 G4int i0;
221 for(i0=0; i0<50; i0++)
222 {
223 it = i0;
224 // G4cout << " SelectExitChannel " << it << " " << random << " " << running[i0]/sum << " " << running[i0] << G4endl; //GDEB
225 if(random < running[i0]/sum) break;
226 }
227 }
228//debug: it = 1;
229// G4cout << " SelectExitChannel " << it << " " << sum << G4endl; //GDEB
230 return it;
231}
virtual G4ParticleHPVector * GetXsec()

References G4UniformRand, G4ParticleHPInelasticCompFS::GetXsec(), G4INCL::Math::max(), and G4ParticleHPInelasticCompFS::theXsection.

Referenced by G4ParticleHPInelasticCompFS::CompositeApply().

◆ 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

◆ 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().

◆ two_body_reaction()

void G4ParticleHPInelasticCompFS::two_body_reaction ( G4ReactionProduct proj,
G4ReactionProduct targ,
G4ReactionProduct product,
G4double  resExcitationEnergy 
)
privateinherited

Definition at line 816 of file G4ParticleHPInelasticCompFS.cc.

820{
821 //CMS system:
822 G4ReactionProduct theCMS= *proj+ *targ;
823
824 //Residual definition:
825 G4int resZ=(G4int)(proj->GetDefinition()->GetPDGCharge()+targ->GetDefinition()->GetPDGCharge()-product->GetDefinition()->GetPDGCharge()+0.1);
826 G4int resA=proj->GetDefinition()->GetBaryonNumber()+targ->GetDefinition()->GetBaryonNumber()-product->GetDefinition()->GetBaryonNumber();
827 G4ReactionProduct theResidual;
828 theResidual.SetDefinition(G4IonTable::GetIonTable()->GetIon(resZ,resA,0.0));
829
830 //CMS system:
831 G4ReactionProduct theCMSproj;
832 G4ReactionProduct theCMStarg;
833 theCMSproj.Lorentz(*proj,theCMS);
834 theCMStarg.Lorentz(*targ,theCMS);
835 //final Momentum in the CMS:
836 G4double totE=std::sqrt(theCMSproj.GetMass()*theCMSproj.GetMass()+theCMSproj.GetTotalMomentum()*theCMSproj.GetTotalMomentum())+std::sqrt(theCMStarg.GetMass()*theCMStarg.GetMass()+theCMStarg.GetTotalMomentum()*theCMStarg.GetTotalMomentum());
837 G4double prodmass=product->GetMass();
838 G4double resmass=theResidual.GetMass()+resExcitationEnergy;
839 G4double fmomsquared=1./4./totE/totE*(totE*totE-(prodmass-resmass)*(prodmass-resmass))*(totE*totE-(prodmass+resmass)*(prodmass+resmass));
840 G4double fmom=0;
841 if(fmomsquared>0){
842 fmom=std::sqrt(fmomsquared);
843 }
844
845 //random (isotropic direction):
846 G4double cosTh = 2.*G4UniformRand()-1.;
848 G4double theta = std::acos(cosTh);
849 G4double sinth = std::sin(theta);
850 product->SetMomentum(fmom*sinth*std::cos(phi),fmom*sinth*std::sin(phi),fmom*cosTh); //CMS
851 product->SetTotalEnergy(std::sqrt(prodmass*prodmass+fmom*fmom)); //CMS
852 //Back to the LAB system:
853 product->Lorentz(*product,-1.*theCMS);
854
855}
void SetTotalEnergy(const G4double en)
G4double GetTotalMomentum() const

References G4UniformRand, G4ParticleDefinition::GetBaryonNumber(), G4ReactionProduct::GetDefinition(), G4IonTable::GetIon(), G4IonTable::GetIonTable(), G4ReactionProduct::GetMass(), G4ParticleDefinition::GetPDGCharge(), G4ReactionProduct::GetTotalMomentum(), G4ReactionProduct::Lorentz(), G4ReactionProduct::SetDefinition(), G4ReactionProduct::SetMomentum(), G4ReactionProduct::SetTotalEnergy(), and CLHEP::twopi.

Referenced by G4ParticleHPInelasticCompFS::CompositeApply().

◆ use_nresp71_model()

G4bool G4ParticleHPInelasticCompFS::use_nresp71_model ( const G4ParticleDefinition aDefinition,
const G4int  it,
const G4ReactionProduct theTarget,
G4ReactionProduct boosted 
)
privateinherited

Definition at line 858 of file G4ParticleHPInelasticCompFS.cc.

859{
860 if ( aDefinition == G4Neutron::Definition() ) // If the outgoing particle is a neutron...
861 {
862 // LR: flag LR in ENDF. It indicates whether there is breakup of the residual nucleus or not.
863 // it: exit channel (index of the carbon excited state)
864
865 //if ( (G4int)(theBaseZ+0.1) == 6 ) G4cout << "LR[" << it << "] = " << LR[it] << G4endl;
866
867 // Added by A. R. GarcĂ­a (CIEMAT) to include the physics of C(N,N'3A) reactions from NRESP71.
868
869 if ( LR[it] > 0 ) // If there is breakup of the residual nucleus LR(flag LR in ENDF)>0 (i.e. Z=6 MT=52-91 (it=MT-50)).
870 {
871 // Defining carbon as the target in the reference frame at rest.
872 G4ReactionProduct theCarbon(theTarget);
873
874 theCarbon.SetMomentum(G4ThreeVector());
875 theCarbon.SetKineticEnergy(0.);
876
877 // Creating four reaction products.
878 G4ReactionProduct theProds[4];
879
880 // Applying C(N,N'3A) reaction mechanisms in the target rest frame.
881 if ( it == 41 )
882 {
883 // QI=QM=-7.275 MeV for C-0(N,N')C-C(3A) in ENDF/B-VII.1. This is not the value of the QI of the first step according to the model. So we don't take it. Instead, we set the one we have calculated: QI=(mn+m12C)-(ma+m9Be+Ex9Be)=-8.130 MeV.
884 nresp71_model.ApplyMechanismI_NBeA2A(boosted, theCarbon, theProds, -8.130/*QI[it]*/); // N+C --> A[0]+9BE* | 9BE* --> N[1]+8BE | 8BE --> 2*A[2,3].
885 //printf("- QI=%f\n", QI[it]);
886 }
887 else
888 {
889 nresp71_model.ApplyMechanismII_ACN2A(boosted, theCarbon, theProds, QI[it]); // N+C --> N'[0]+C* | C* --> A[1]+8BE | 8BE --> 2*A[2,3].
890 }
891
892 //printf("it=%d qi=%f \n", it, QI[it]);
893
894 // Returning to the reference frame where the target was in motion.
895 for ( G4int j=0; j<4; j++ )
896 {
897 theProds[j].Lorentz(theProds[j], -1.*theTarget);
898 theResult.Get()->AddSecondary(new G4DynamicParticle(theProds[j].GetDefinition(), theProds[j].GetMomentum()), secID);
899 }
900
901 /*G4double EN0 = theNeutron.GetKineticEnergy();
902 G4double EN1 = theProds[0].GetKineticEnergy();
903
904 G4double EA1 = theProds[1].GetKineticEnergy();
905 G4double EA2 = theProds[2].GetKineticEnergy();
906 G4double EA3 = theProds[3].GetKineticEnergy();
907
908 printf("Q=%f\n", EN1+EA1+EA2+EA3-EN0);*/
909
910 // Killing the primary neutron.
912
913 return true;
914 }
915 }
916 else if ( aDefinition == G4Alpha::Definition() ) // If the outgoing particle is an alpha, ...
917 {
918 // Added by A. R. GarcĂ­a (CIEMAT) to include the physics of C(N,A)9BE reactions from NRESP71.
919
920 if ( LR[it] == 0 ) // If Z=6, an alpha particle is emitted and there is no breakup of the residual nucleus LR(flag LR in ENDF)==0.
921 {
922 // Defining carbon as the target in the reference frame at rest.
923 G4ReactionProduct theCarbon(theTarget);
924
925 theCarbon.SetMomentum(G4ThreeVector());
926 theCarbon.SetKineticEnergy(0.);
927
928 // Creating four reaction products.
929 G4ReactionProduct theProds[2];
930
931 // Applying C(N,A)9BE reaction mechanism.
932 nresp71_model.ApplyMechanismABE(boosted, theCarbon, theProds); // N+C --> A[0]+9BE[1].
933
934 //G4DynamicParticle *theSec;
935 for ( G4int j=0; j<2; j++ )
936 {
937 // Returning to the system of reference where the target was in motion.
938 theProds[j].Lorentz(theProds[j], -1.*theTarget);
939 theResult.Get()->AddSecondary(new G4DynamicParticle(theProds[j].GetDefinition(), theProds[j].GetMomentum()), secID);
940 }
941
942 // Killing the primary neutron.
944
945 return true;
946 }
947 else
948 {
949 G4Exception("G4ParticleHPInelasticCompFS::CompositeApply()", "G4ParticleInelasticCompFS.cc", FatalException, "Alpha production with LR!=0.");
950 }
951 }
952
953 return false;
954}
@ FatalException
CLHEP::Hep3Vector G4ThreeVector
static G4Alpha * Definition()
Definition: G4Alpha.cc:48
G4int ApplyMechanismABE(G4ReactionProduct &neut, G4ReactionProduct &carb, G4ReactionProduct *theProds)
G4int ApplyMechanismI_NBeA2A(G4ReactionProduct &neut, G4ReactionProduct &carb, G4ReactionProduct *theProds, const G4double QI)
G4int ApplyMechanismII_ACN2A(G4ReactionProduct &neut, G4ReactionProduct &carb, G4ReactionProduct *theProds, const G4double QI)

References G4HadFinalState::AddSecondary(), G4NRESP71M03::ApplyMechanismABE(), G4NRESP71M03::ApplyMechanismI_NBeA2A(), G4NRESP71M03::ApplyMechanismII_ACN2A(), G4Neutron::Definition(), G4Alpha::Definition(), FatalException, G4Exception(), G4Cache< VALTYPE >::Get(), G4ReactionProduct::Lorentz(), G4ParticleHPInelasticCompFS::LR, G4ParticleHPInelasticCompFS::nresp71_model, G4ParticleHPInelasticCompFS::QI, G4ParticleHPFinalState::secID, G4ReactionProduct::SetKineticEnergy(), G4ReactionProduct::SetMomentum(), G4HadFinalState::SetStatusChange(), stopAndKill, and G4ParticleHPFinalState::theResult.

Referenced by G4ParticleHPInelasticCompFS::CompositeApply().

Field Documentation

◆ gammaPath

G4String G4ParticleHPInelasticCompFS::gammaPath
protectedinherited

◆ hasAnyData

G4bool G4ParticleHPFinalState::hasAnyData
protectedinherited

◆ hasFSData

G4bool G4ParticleHPFinalState::hasFSData
protectedinherited

◆ hasXsec

G4bool G4ParticleHPFinalState::hasXsec
protectedinherited

◆ LR

std::vector<G4int> G4ParticleHPInelasticCompFS::LR
protectedinherited

◆ nresp71_model

G4NRESP71M03 G4ParticleHPInelasticCompFS::nresp71_model
privateinherited

◆ QI

std::vector<G4double> G4ParticleHPInelasticCompFS::QI
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(), 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().

◆ theAngularDistribution

G4ParticleHPAngular* G4ParticleHPInelasticCompFS::theAngularDistribution[51]
protectedinherited

◆ theBaseA

G4double G4ParticleHPFinalState::theBaseA
protectedinherited

◆ theBaseM

G4int G4ParticleHPFinalState::theBaseM
protectedinherited

◆ theBaseZ

G4double G4ParticleHPFinalState::theBaseZ
protectedinherited

◆ theEnergyAngData

G4ParticleHPEnAngCorrelation* G4ParticleHPInelasticCompFS::theEnergyAngData[51]
protectedinherited

◆ theEnergyDistribution

G4ParticleHPEnergyDistribution* G4ParticleHPInelasticCompFS::theEnergyDistribution[51]
protectedinherited

◆ theFinalStatePhotons

G4ParticleHPPhotonDist* G4ParticleHPInelasticCompFS::theFinalStatePhotons[51]
protectedinherited

◆ theGammas

G4ParticleHPDeExGammas G4ParticleHPInelasticCompFS::theGammas
protectedinherited

◆ theNames

G4ParticleHPNames G4ParticleHPFinalState::theNames
protectedinherited

◆ theNDLDataA

G4int G4ParticleHPFinalState::theNDLDataA
protectedinherited

◆ theNDLDataM

G4int G4ParticleHPFinalState::theNDLDataM
protectedinherited

◆ theNDLDataZ

G4int G4ParticleHPFinalState::theNDLDataZ
protectedinherited

◆ theProjectile

G4ParticleDefinition* G4ParticleHPFinalState::theProjectile
protectedinherited

◆ 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(), 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().

◆ theXsection

G4ParticleHPVector* G4ParticleHPInelasticCompFS::theXsection[51]
protectedinherited

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