Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4BinaryLightIonReaction Class Reference

#include <G4BinaryLightIonReaction.hh>

+ Inheritance diagram for G4BinaryLightIonReaction:

Public Member Functions

 G4BinaryLightIonReaction (G4VPreCompoundModel *ptr=0)
 
virtual ~G4BinaryLightIonReaction ()
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
void SetPrecompound (G4VPreCompoundModel *ptr)
 
void SetDeExcitation (G4ExcitationHandler *ptr)
 
virtual void ModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair< G4double, G4doubleGetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 
 G4HadronicInteraction (const G4HadronicInteraction &right)=delete
 
const G4HadronicInteractionoperator= (const G4HadronicInteraction &right)=delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 34 of file G4BinaryLightIonReaction.hh.

Constructor & Destructor Documentation

◆ G4BinaryLightIonReaction()

G4BinaryLightIonReaction::G4BinaryLightIonReaction ( G4VPreCompoundModel ptr = 0)

Definition at line 57 of file G4BinaryLightIonReaction.cc.

58: G4HadronicInteraction("Binary Light Ion Cascade"),
59 theProjectileFragmentation(ptr),
60 pA(0),pZ(0), tA(0),tZ(0),spectatorA(0),spectatorZ(0),
61 projectile3dNucleus(0),target3dNucleus(0)
62{
63 if(!ptr) {
66 G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(p);
67 if(!pre) { pre = new G4PreCompoundModel(); }
68 theProjectileFragmentation = pre;
69 }
70 theModel = new G4BinaryCascade(theProjectileFragmentation);
71 theHandler = theProjectileFragmentation->GetExcitationHandler();
72 if ( theBLIR_ID == -1 ) {
73#ifdef G4MULTITHREADED
74 G4MUTEXLOCK(&G4BinaryLightIonReaction::BLIRMutex);
75 if ( theBLIR_ID == -1 ) {
76#endif
77 theBLIR_ID = G4PhysicsModelCatalog::Register("Binary Light Ion Reaction");
78#ifdef G4MULTITHREADED
79 }
80 G4MUTEXUNLOCK(&G4BinaryLightIonReaction::BLIRMutex);
81#endif
82 }
83
84
85 debug_G4BinaryLightIonReactionResults=std::getenv("debug_G4BinaryLightIonReactionResults")!=0;
86}
#define G4MUTEXLOCK(mutex)
Definition: G4Threading.hh:251
#define G4MUTEXUNLOCK(mutex)
Definition: G4Threading.hh:254
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
static G4int Register(const G4String &)
G4ExcitationHandler * GetExcitationHandler() const

◆ ~G4BinaryLightIonReaction()

G4BinaryLightIonReaction::~G4BinaryLightIonReaction ( )
virtual

Definition at line 88 of file G4BinaryLightIonReaction.cc.

89{}

Member Function Documentation

◆ ApplyYourself()

G4HadFinalState * G4BinaryLightIonReaction::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus theNucleus 
)
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 106 of file G4BinaryLightIonReaction.cc.

108{
109 if(std::getenv("BLICDEBUG") ) G4cerr << " ######### Binary Light Ion Reaction starts ######### " << G4endl;
110 G4ping debug("debug_G4BinaryLightIonReaction");
111 pA=aTrack.GetDefinition()->GetBaryonNumber();
112 pZ=G4lrint(aTrack.GetDefinition()->GetPDGCharge()/eplus);
113 tA=targetNucleus.GetA_asInt();
114 tZ=targetNucleus.GetZ_asInt();
115 G4double timePrimary = aTrack.GetGlobalTime();
116 G4LorentzVector mom(aTrack.Get4Momentum());
117 //G4cout << "proj mom : " << mom << G4endl;
118 G4LorentzRotation toBreit(mom.boostVector());
119
120 G4bool swapped=SetLighterAsProjectile(mom, toBreit);
121 //G4cout << "after swap, swapped? / mom " << swapped << " / " << mom <<G4endl;
122 G4ReactionProductVector * result = 0;
123 G4ReactionProductVector * cascaders=0; //new G4ReactionProductVector;
124// G4double m_nucl(0); // to check energy balance
125
126
127 // G4double m1=G4ParticleTable::GetParticleTable()->GetIonTable()->GetIonMass(pZ,pA);
128 // G4cout << "Entering the decision point "
129 // << (mom.t()-mom.mag())/pA << " "
130 // << pA<<" "<< pZ<<" "
131 // << tA<<" "<< tZ<<G4endl
132 // << " "<<mom.t()-mom.mag()<<" "
133 // << mom.t()- m1<<G4endl;
134 if( (mom.t()-mom.mag())/pA < 50*MeV )
135 {
136 // G4cout << "Using pre-compound only, E= "<<mom.t()-mom.mag()<<G4endl;
137 // m_nucl = mom.mag();
138 cascaders=FuseNucleiAndPrompound(mom);
139 if( !cascaders )
140 {
141
142 // abort!! happens for too low energy for nuclei to fuse
143
144 theResult.Clear();
145 theResult.SetStatusChange(isAlive);
146 theResult.SetEnergyChange(aTrack.GetKineticEnergy());
147 theResult.SetMomentumChange(aTrack.Get4Momentum().vect().unit());
148 return &theResult;
149 }
150 }
151 else
152 {
153 result=Interact(mom,toBreit);
154
155 if(! result )
156 {
157 // abort!!
158
159 G4cerr << "G4BinaryLightIonReaction no final state for: " << G4endl;
160 G4cerr << " Primary " << aTrack.GetDefinition()
161 << ", (A,Z)=(" << aTrack.GetDefinition()->GetBaryonNumber()
162 << "," << aTrack.GetDefinition()->GetPDGCharge()/eplus << ") "
163 << ", kinetic energy " << aTrack.GetKineticEnergy()
164 << G4endl;
165 G4cerr << " Target nucleus (A,Z)=("
166 << (swapped?pA:tA) << ","
167 << (swapped?pZ:tZ) << ")" << G4endl;
168 G4cerr << " if frequent, please submit above information as bug report"
169 << G4endl << G4endl;
170
171 theResult.Clear();
172 theResult.SetStatusChange(isAlive);
173 theResult.SetEnergyChange(aTrack.GetKineticEnergy());
174 theResult.SetMomentumChange(aTrack.Get4Momentum().vect().unit());
175 return &theResult;
176 }
177
178 // Calculate excitation energy,
179 G4double theStatisticalExEnergy = GetProjectileExcitation();
180
181
182 pInitialState = mom;
183 //G4cout << "BLIC: pInitialState from aTrack : " << pInitialState;
184 pInitialState.setT(pInitialState.getT() +
186 //G4cout << "BLIC: target nucleus added : " << pInitialState << G4endl;
187
188 delete target3dNucleus;target3dNucleus=0;
189 delete projectile3dNucleus;projectile3dNucleus=0;
190
192
193 cascaders = new G4ReactionProductVector;
194
195 G4LorentzVector pspectators=SortResult(result,spectators,cascaders);
196 // this also sets spectatorA and spectatorZ
197
198 // pFinalState=std::accumulate(cascaders->begin(),cascaders->end(),pFinalState,ReactionProduct4Mom);
199
200 std::vector<G4ReactionProduct *>::iterator iter;
201
202 // G4cout << "pInitialState, pFinalState / pspectators"<< pInitialState << " / " << pFinalState << " / " << pspectators << G4endl;
203 // if ( spectA-spectatorA !=0 || spectZ-spectatorZ !=0)
204 // {
205 // G4cout << "spect Nucl != spectators: nucl a,z; spect a,z" <<
206 // spectatorA <<" "<< spectatorZ <<" ; " << spectA <<" "<< spectZ << G4endl;
207 // }
208 delete result;
209 result=0;
210 G4LorentzVector momentum(pInitialState-pFinalState);
211 G4int loopcount(0);
212 //G4cout << "BLIC: momentum, pspectators : " << momentum << " / " << pspectators << G4endl;
213 while (std::abs(momentum.e()-pspectators.e()) > 10*MeV) /* Loop checking, 31.08.2015, G.Folger */
214 // see if on loopcount
215 {
216 G4LorentzVector pCorrect(pInitialState-pspectators);
217 //G4cout << "BLIC:: BIC nonconservation? (pInitialState-pFinalState) / spectators :" << momentum << " / " << pspectators << "pCorrect "<< pCorrect<< G4endl;
218 // Correct outgoing casacde particles.... to have momentum of (initial state - spectators)
219 G4bool EnergyIsCorrect=EnergyAndMomentumCorrector(cascaders, pCorrect);
220 if ( ! EnergyIsCorrect && debug_G4BinaryLightIonReactionResults)
221 {
222 G4cout << "Warning - G4BinaryLightIonReaction E/P correction for cascaders failed" << G4endl;
223 }
224 pFinalState=G4LorentzVector(0,0,0,0);
225 for(iter=cascaders->begin(); iter!=cascaders->end(); iter++)
226 {
227 pFinalState += G4LorentzVector( (*iter)->GetMomentum(), (*iter)->GetTotalEnergy() );
228 }
229 momentum=pInitialState-pFinalState;
230 if (++loopcount > 10 )
231 {
232 if ( momentum.vect().mag() - momentum.e()> 10*keV )
233 {
234 G4cerr << "G4BinaryLightIonReaction.cc: Cannot correct 4-momentum of cascade particles" << G4endl;
235 throw G4HadronicException(__FILE__, __LINE__, "G4BinaryCasacde::ApplyCollision()");
236 } else {
237 break;
238 }
239
240 }
241 }
242
243 if (spectatorA > 0 )
244 {
245 // check spectator momentum
246 if ( momentum.vect().mag() - momentum.e()< 10*keV )
247 {
248 // DeExciteSpectatorNucleus() also handles also case of A=1, Z=0,1
249 DeExciteSpectatorNucleus(spectators, cascaders, theStatisticalExEnergy, momentum);
250
251 } else { // momentum non-conservation --> fail
252 for (iter=spectators->begin();iter!=spectators->end();iter++)
253 {
254 delete *iter;
255 }
256 delete spectators;
257 for(iter=cascaders->begin(); iter!=cascaders->end(); iter++)
258 {
259 delete *iter;
260 }
261 delete cascaders;
262
263 G4cout << "G4BinaryLightIonReaction.cc: mom check: " << momentum
264 << " 3.mag "<< momentum.vect().mag() << G4endl
265 << " .. pInitialState/pFinalState/spectators " << pInitialState <<" "
266 << pFinalState << " " << pspectators << G4endl
267 << " .. A,Z " << spectatorA <<" "<< spectatorZ << G4endl;
268 G4cout << "G4BinaryLightIonReaction invalid final state for: " << G4endl;
269 G4cout << " Primary " << aTrack.GetDefinition()
270 << ", (A,Z)=(" << aTrack.GetDefinition()->GetBaryonNumber()
271 << "," << aTrack.GetDefinition()->GetPDGCharge()/eplus << ") "
272 << ", kinetic energy " << aTrack.GetKineticEnergy()
273 << G4endl;
274 G4cout << " Target nucleus (A,Z)=(" << targetNucleus.GetA_asInt()
275 << "," << targetNucleus.GetZ_asInt() << ")" << G4endl;
276 G4cout << " if frequent, please submit above information as bug report"
277 << G4endl << G4endl;
278#ifdef debug_G4BinaryLightIonReaction
280 ed << "G4BinaryLightIonreaction: Terminate for above error" << G4endl;
281 G4Exception("G4BinaryLightIonreaction::ApplyYourSelf()", "BLIC001", FatalException,
282 ed);
283
284#endif
285 theResult.Clear();
286 theResult.SetStatusChange(isAlive);
287 theResult.SetEnergyChange(aTrack.GetKineticEnergy());
288 theResult.SetMomentumChange(aTrack.Get4Momentum().vect().unit());
289 return &theResult;
290 }
291 } else { // no spectators
292 delete spectators;
293 }
294 }
295 // Rotate to lab
297 toZ.rotateZ(-1*mom.phi());
298 toZ.rotateY(-1*mom.theta());
299 G4LorentzRotation toLab(toZ.inverse());
300
301 // Fill the particle change, while rotating. Boost from projectile breit-frame in case we swapped.
302 // theResult.Clear();
303 theResult.Clear();
304 theResult.SetStatusChange(stopAndKill);
305 G4LorentzVector ptot(0);
306 G4ReactionProductVector::iterator iter;
307 #ifdef debug_BLIR_result
308 G4LorentzVector p_raw;
309 #endif
310 //G4int i=0;
311
312 for(iter=cascaders->begin(); iter!=cascaders->end(); iter++)
313 {
314 if((*iter)->GetNewlyAdded())
315 {
316 G4DynamicParticle * aNewDP =
317 new G4DynamicParticle((*iter)->GetDefinition(),
318 (*iter)->GetTotalEnergy(),
319 (*iter)->GetMomentum() );
320 G4LorentzVector tmp = aNewDP->Get4Momentum();
321 #ifdef debug_BLIR_result
322 p_raw+= tmp;
323 #endif
324 if(swapped)
325 {
326 tmp*=toBreit.inverse();
327 tmp.setVect(-tmp.vect());
328 }
329 tmp *= toLab;
330 aNewDP->Set4Momentum(tmp);
331 G4HadSecondary aNew = G4HadSecondary(aNewDP);
332 G4double time = 0; //(*iter)->GetCreationTime();
333 //if(time < 0.0) { time = 0.0; }
334 aNew.SetTime(timePrimary + time);
335 aNew.SetCreatorModelType((*iter)->GetCreatorModel());
336
337 theResult.AddSecondary(aNew);
338 ptot += tmp;
339 //G4cout << "BLIC: Secondary " << aNew->GetDefinition()->GetParticleName()
340 // <<" "<< aNew->GetMomentum()<<" "<< aNew->GetTotalEnergy() << G4endl;
341 }
342 delete *iter;
343 }
344 delete cascaders;
345
346#ifdef debug_BLIR_result
347 //G4cout << "Result analysis, secondaries " << theResult.GetNumberOfSecondaries() << G4endl;
348 //G4cout << "p_tot_raw " << p_raw << " sum p final " << ptot << G4endl;
350 GetIonMass(targetNucleus.GetZ_asInt(),targetNucleus.GetA_asInt());
351 // delete? tZ=targetNucleus.GetZ_asInt();
352
353 //G4cout << "BLIC Energy conservation initial/primary/nucleus/final/delta(init-final) "
354 // << aTrack.GetTotalEnergy() + m_nucl <<" "<< aTrack.GetTotalEnergy() <<" "<< m_nucl <<" "<<ptot.e()
355 // <<" "<< aTrack.GetTotalEnergy() + m_nucl - ptot.e() << G4endl;
356 G4cout << "BLIC momentum conservation " << aTrack.Get4Momentum()+ G4LorentzVector(m_nucl)
357 << " ptot " << ptot << " delta " << aTrack.Get4Momentum()+ G4LorentzVector(m_nucl) - ptot
358 << " 3mom.mag() " << (aTrack.Get4Momentum()+ G4LorentzVector(m_nucl) - ptot).vect().mag() << G4endl;
359#endif
360
361 if(std::getenv("BLICDEBUG") ) G4cerr << " ######### Binary Light Ion Reaction number ends ######### " << G4endl;
362
363 return &theResult;
364}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
@ isAlive
@ stopAndKill
CLHEP::HepLorentzVector G4LorentzVector
std::vector< G4ReactionProduct * > G4ReactionProductVector
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
Hep3Vector unit() const
HepLorentzRotation & rotateY(double delta)
HepLorentzRotation & rotateZ(double delta)
HepLorentzRotation inverse() const
Hep3Vector vect() const
void setVect(const Hep3Vector &)
G4LorentzVector Get4Momentum() const
void Set4Momentum(const G4LorentzVector &momentum)
void SetStatusChange(G4HadFinalStateStatus aS)
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
void SetEnergyChange(G4double anEnergy)
void SetMomentumChange(const G4ThreeVector &aV)
const G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
G4double GetGlobalTime() const
void SetTime(G4double aT)
void SetCreatorModelType(G4int idx)
G4double GetIonMass(G4int Z, G4int A, G4int L=0, G4int lvl=0) const
Definition: G4IonTable.cc:1517
G4double GetPDGCharge() const
G4IonTable * GetIonTable() const
static G4ParticleTable * GetParticleTable()
Definition: G4ping.hh:35
int G4lrint(double ad)
Definition: templates.hh:134

◆ ModelDescription()

void G4BinaryLightIonReaction::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 91 of file G4BinaryLightIonReaction.cc.

92{
93 outFile << "G4Binary Light Ion Cascade is an intra-nuclear cascade model\n"
94 << "using G4BinaryCasacde to model the interaction of a light\n"
95 << "nucleus with a nucleus.\n"
96 << "The lighter of the two nuclei is treated like a set of projectiles\n"
97 << "which are transported simultaneously through the heavier nucleus.\n";
98}

◆ SetDeExcitation()

void G4BinaryLightIonReaction::SetDeExcitation ( G4ExcitationHandler ptr)
inline

Definition at line 78 of file G4BinaryLightIonReaction.hh.

79{
80 theProjectileFragmentation->SetExcitationHandler(ptr);
81 theHandler = ptr;
82}
void SetExcitationHandler(G4ExcitationHandler *ptr)

◆ SetPrecompound()

void G4BinaryLightIonReaction::SetPrecompound ( G4VPreCompoundModel ptr)
inline

Definition at line 73 of file G4BinaryLightIonReaction.hh.

74{
75 if(ptr) { theProjectileFragmentation = ptr; }
76 theHandler = theProjectileFragmentation->GetExcitationHandler();
77}

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