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

#include <G4ParticleHPAngular.hh>

Data Structures

struct  toBeCached
 

Public Member Functions

 G4ParticleHPAngular ()
 
G4double GetTargetMass ()
 
void Init (std::istream &aDataFile)
 
void SampleAndUpdate (G4ReactionProduct &anIncidentParticle)
 
void SetProjectileRP (const G4ReactionProduct &anIncidentParticleRP)
 
void SetTarget (const G4ReactionProduct &aTarget)
 
 ~G4ParticleHPAngular ()
 

Private Attributes

G4Cache< toBeCachedfCache
 
G4int frameFlag
 
G4double targetMass
 
G4int theAngularDistributionType
 
G4ParticleHPLegendreStoretheCoefficients
 
G4bool theIsoFlag
 
G4ParticleHPPartialtheProbArray
 

Detailed Description

Definition at line 45 of file G4ParticleHPAngular.hh.

Constructor & Destructor Documentation

◆ G4ParticleHPAngular()

G4ParticleHPAngular::G4ParticleHPAngular ( )
inline

Definition at line 56 of file G4ParticleHPAngular.hh.

57 {
58 theIsoFlag = true;
60 theProbArray = 0;
61 toBeCached val;
62 fCache.Put( val );
63
65 frameFlag = 0;
66 targetMass = 0.0;
67 }
G4ParticleHPPartial * theProbArray
G4Cache< toBeCached > fCache
G4ParticleHPLegendreStore * theCoefficients

References fCache, frameFlag, targetMass, theAngularDistributionType, theCoefficients, theIsoFlag, and theProbArray.

◆ ~G4ParticleHPAngular()

G4ParticleHPAngular::~G4ParticleHPAngular ( )
inline

Definition at line 69 of file G4ParticleHPAngular.hh.

70 {
71 delete theCoefficients;
72 delete theProbArray;
73 }

References theCoefficients, and theProbArray.

Member Function Documentation

◆ GetTargetMass()

G4double G4ParticleHPAngular::GetTargetMass ( )
inline

Definition at line 89 of file G4ParticleHPAngular.hh.

90 {
91 return targetMass;
92 }

References targetMass.

Referenced by G4ParticleHPInelasticBaseFS::BaseApply().

◆ Init()

void G4ParticleHPAngular::Init ( std::istream &  aDataFile)

Definition at line 41 of file G4ParticleHPAngular.cc.

42{
43// G4cout << "here we are entering the Angular Init"<<G4endl;
45 aDataFile >> frameFlag;
47 {
48 theIsoFlag = true;
49 }
51 {
52 theIsoFlag = false;
53 G4int nEnergy;
54 aDataFile >> nEnergy;
57 G4double temp, energy;
58 G4int tempdep, nLegendre;
59 G4int i, ii;
60 for (i=0; i<nEnergy; i++)
61 {
62 aDataFile >> temp >> energy >> tempdep >> nLegendre;
63 energy *=eV;
64 theCoefficients->Init(i, energy, nLegendre);
66 G4double coeff=0;
67 for(ii=0; ii<nLegendre; ii++)
68 {
69 aDataFile >> coeff;
70 theCoefficients->SetCoeff(i, ii+1, coeff);
71 }
72 }
73 }
75 {
76 theIsoFlag = false;
77 G4int nEnergy;
78 aDataFile >> nEnergy;
79 theProbArray = new G4ParticleHPPartial(nEnergy, nEnergy);
81 G4double temp, energy;
82 G4int tempdep;
83 for(G4int i=0; i<nEnergy; i++)
84 {
85 aDataFile >> temp >> energy >> tempdep;
86 energy *= eV;
87 theProbArray->SetT(i, temp);
89 theProbArray->InitData(i, aDataFile);
90 }
91 }
92 else
93 {
94 theIsoFlag = false;
95 G4cout << "unknown distribution found for Angular: "<< theAngularDistributionType << G4endl;
96 throw G4HadronicException(__FILE__, __LINE__, "unknown distribution needs implementation!!!");
97 }
98}
static constexpr double eV
Definition: G4SIunits.hh:201
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
void SetCoeff(G4int i, G4int l, G4double coeff)
void Init(G4int i, G4double e, G4int n)
void SetTemperature(G4int i, G4double temp)
void InitInterpolation(std::istream &aDataFile)
void SetT(G4int i, G4double x)
void SetX(G4int i, G4double x)
void InitData(G4int i, std::istream &aDataFile, G4double unit=1.)
void InitInterpolation(G4int i, std::istream &aDataFile)
G4double energy(const ThreeVector &p, const G4double m)

References G4INCL::KinematicsUtils::energy(), eV, frameFlag, G4cout, G4endl, G4ParticleHPLegendreStore::Init(), G4ParticleHPPartial::InitData(), G4ParticleHPPartial::InitInterpolation(), G4ParticleHPLegendreStore::InitInterpolation(), G4ParticleHPLegendreStore::SetCoeff(), G4ParticleHPPartial::SetT(), G4ParticleHPLegendreStore::SetTemperature(), G4ParticleHPPartial::SetX(), targetMass, theAngularDistributionType, theCoefficients, theIsoFlag, and theProbArray.

Referenced by G4FissionLibrary::Init(), G4ParticleHPFSFissionFS::Init(), G4ParticleHPInelasticCompFS::Init(), G4ParticleHPFissionBaseFS::Init(), and G4ParticleHPInelasticBaseFS::Init().

◆ SampleAndUpdate()

void G4ParticleHPAngular::SampleAndUpdate ( G4ReactionProduct anIncidentParticle)

Definition at line 100 of file G4ParticleHPAngular.cc.

101{
102
103 //********************************************************************
104 //EMendoza -> sampling can be isotropic in LAB or in CMS
105 /*
106 if(theIsoFlag)
107 {
108// G4cout << "Angular result "<<aHadron.GetTotalMomentum()<<" ";
109// @@@ add code for isotropic emission in CMS.
110 G4double costheta = 2.*G4UniformRand()-1;
111 G4double theta = std::acos(costheta);
112 G4double phi = twopi*G4UniformRand();
113 G4double sinth = std::sin(theta);
114 G4double en = aHadron.GetTotalMomentum();
115 G4ThreeVector temp(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
116 aHadron.SetMomentum( temp );
117 aHadron.Lorentz(aHadron, -1.*theTarget);
118 }
119 else
120 {
121 */
122 //********************************************************************
123 if(frameFlag == 1) // LAB
124 {
125 G4double en = aHadron.GetTotalMomentum();
126 G4ReactionProduct boosted;
127 boosted.Lorentz(*fCache.Get().theProjectileRP, *fCache.Get().theTarget);
128 G4double kineticEnergy = boosted.GetKineticEnergy();
129 G4double cosTh = 0.0;
130 //********************************************************************
131 //EMendoza --> sampling can be also isotropic
132 /*
133 if(theAngularDistributionType == 1) cosTh = theCoefficients->SampleMax(kineticEnergy);
134 if(theAngularDistributionType == 2) cosTh = theProbArray->Sample(kineticEnergy);
135 */
136 //********************************************************************
137 if(theIsoFlag){cosTh =2.*G4UniformRand()-1;}
138 else if(theAngularDistributionType == 1) {cosTh = theCoefficients->SampleMax(kineticEnergy);}
139 else if(theAngularDistributionType == 2) {cosTh = theProbArray->Sample(kineticEnergy);}
140 else{
141 G4cout << "unknown distribution found for Angular: "<< theAngularDistributionType << G4endl;
142 throw G4HadronicException(__FILE__, __LINE__, "unknown distribution needs implementation!!!");
143 }
144 //********************************************************************
145 G4double theta = std::acos(cosTh);
147 G4double sinth = std::sin(theta);
148 G4ThreeVector temp(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
149 aHadron.SetMomentum( temp );
150 }
151 else if(frameFlag == 2) // costh in CMS
152 {
153 G4ReactionProduct boostedN;
154 boostedN.Lorentz(*fCache.Get().theProjectileRP, *fCache.Get().theTarget);
155 G4double kineticEnergy = boostedN.GetKineticEnergy();
156
157 G4double cosTh = 0.0;
158 //********************************************************************
159 //EMendoza --> sampling can be also isotropic
160 /*
161 if(theAngularDistributionType == 1) cosTh = theCoefficients->SampleMax(kineticEnergy);
162 if(theAngularDistributionType == 2) cosTh = theProbArray->Sample(kineticEnergy);
163 */
164 //********************************************************************
165 if(theIsoFlag){cosTh =2.*G4UniformRand()-1;}
166 else if(theAngularDistributionType == 1) {cosTh = theCoefficients->SampleMax(kineticEnergy);}
167 else if(theAngularDistributionType == 2) {cosTh = theProbArray->Sample(kineticEnergy);}
168 else{
169 G4cout << "unknown distribution found for Angular: "<< theAngularDistributionType << G4endl;
170 throw G4HadronicException(__FILE__, __LINE__, "unknown distribution needs implementation!!!");
171 }
172 //********************************************************************
173//080612TK bug fix contribution from Benoit Pirard and Laurent Desorgher (Univ. Bern)
174/*
175 if(theAngularDistributionType == 1) // LAB
176 {
177 G4double en = aHadron.GetTotalMomentum();
178 G4ReactionProduct boosted;
179 boosted.Lorentz(theProjectile, theTarget);
180 G4double kineticEnergy = boosted.GetKineticEnergy();
181 G4double cosTh = theCoefficients->SampleMax(kineticEnergy);
182 G4double theta = std::acos(cosTh);
183 G4double phi = twopi*G4UniformRand();
184 G4double sinth = std::sin(theta);
185 G4ThreeVector temp(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
186 aHadron.SetMomentum( temp );
187 }
188 else if(theAngularDistributionType == 2) // costh in CMS {
189 }
190*/
191
192// G4ReactionProduct boostedN;
193// boostedN.Lorentz(theProjectile, theTarget);
194// G4double kineticEnergy = boostedN.GetKineticEnergy();
195// G4double cosTh = theProbArray->Sample(kineticEnergy);
196
197 G4double theta = std::acos(cosTh);
199 G4double sinth = std::sin(theta);
200
201 G4ThreeVector temp(sinth*std::cos(phi), sinth*std::sin(phi), std::cos(theta) ); //CMS
202
203//080612TK bug fix contribution from Benoit Pirard and Laurent Desorgher (Univ. Bern) #5
204/*
205 G4double en = aHadron.GetTotalEnergy(); // Target rest
206
207 // get trafo from Target rest frame to CMS
208 G4ReactionProduct boostedT;
209 boostedT.Lorentz(theTarget, theTarget);
210
211 G4ThreeVector the3IncidentParticle = boostedN.GetMomentum();
212 G4double nEnergy = boostedN.GetTotalEnergy();
213 G4ThreeVector the3Target = boostedT.GetMomentum();
214 G4double tEnergy = boostedT.GetTotalEnergy();
215 G4double totE = nEnergy+tEnergy;
216 G4ThreeVector the3trafo = -the3Target-the3IncidentParticle;
217 G4ReactionProduct trafo; // for transformation from CMS to target rest frame
218 trafo.SetMomentum(the3trafo);
219 G4double cmsMom = std::sqrt(the3trafo*the3trafo);
220 G4double sqrts = std::sqrt((totE-cmsMom)*(totE+cmsMom));
221 trafo.SetMass(sqrts);
222 trafo.SetTotalEnergy(totE);
223
224 G4double gamma = trafo.GetTotalEnergy()/trafo.GetMass();
225 G4double cosalpha = temp*trafo.GetMomentum()/trafo.GetTotalMomentum()/temp.mag();
226 G4double fac = cosalpha*trafo.GetTotalMomentum()/trafo.GetMass();
227 fac*=gamma;
228
229 G4double mom;
230// For G4FPE_DEBUG ON
231 G4double mom2 = ( en*fac*en*fac -
232 (fac*fac - gamma*gamma)*
233 (en*en - gamma*gamma*aHadron.GetMass()*aHadron.GetMass())
234 );
235 if ( mom2 > 0.0 )
236 mom = std::sqrt( mom2 );
237 else
238 mom = 0.0;
239
240 mom = -en*fac - mom;
241 mom /= (fac*fac-gamma*gamma);
242 temp = mom*temp;
243
244 aHadron.SetMomentum( temp ); // now all in CMS
245 aHadron.SetTotalEnergy( std::sqrt( mom*mom + aHadron.GetMass()*aHadron.GetMass() ) );
246 aHadron.Lorentz(aHadron, trafo); // now in target rest frame
247*/
248 // Determination of the hadron kinetic energy in CMS
249 // aHadron.GetKineticEnergy() is actually the residual kinetic energy in CMS (or target frame)
250 // kineticEnergy is incident neutron kinetic energy in CMS (or target frame)
251 G4double QValue = aHadron.GetKineticEnergy() - kineticEnergy;
252 G4double A1 = fCache.Get().theTarget->GetMass()/boostedN.GetMass();
253 G4double A1prim = aHadron.GetMass()/ boostedN.GetMass();
254 G4double kinE = (A1+1-A1prim)/(A1+1)/(A1+1)*(A1*kineticEnergy+(1+A1)*QValue);
255 G4double totalE = kinE + aHadron.GetMass();
256 G4double mom2 = totalE*totalE - aHadron.GetMass()*aHadron.GetMass();
257 G4double mom;
258 if ( mom2 > 0.0 ) mom = std::sqrt( mom2 );
259 else mom = 0.0;
260
261 aHadron.SetMomentum( mom*temp ); // Set momentum in CMS
262 aHadron.SetKineticEnergy(kinE); // Set kinetic energy in CMS
263
264 // get trafo from Target rest frame to CMS
265 G4ReactionProduct boostedT;
266 boostedT.Lorentz(*fCache.Get().theTarget, *fCache.Get().theTarget);
267
268 G4ThreeVector the3IncidentParticle = boostedN.GetMomentum();
269 G4double nEnergy = boostedN.GetTotalEnergy();
270 G4ThreeVector the3Target = boostedT.GetMomentum();
271 G4double tEnergy = boostedT.GetTotalEnergy();
272 G4double totE = nEnergy+tEnergy;
273 G4ThreeVector the3trafo = -the3Target-the3IncidentParticle;
274 G4ReactionProduct trafo; // for transformation from CMS to target rest frame
275 trafo.SetMomentum(the3trafo);
276 G4double cmsMom = std::sqrt(the3trafo*the3trafo);
277 G4double sqrts = std::sqrt((totE-cmsMom)*(totE+cmsMom));
278 trafo.SetMass(sqrts);
279 trafo.SetTotalEnergy(totE);
280
281 aHadron.Lorentz(aHadron, trafo);
282
283 }
284 else
285 {
286 throw G4HadronicException(__FILE__, __LINE__, "Tried to sample non isotropic neutron angular");
287 }
288 aHadron.Lorentz(aHadron, -1.*(*fCache.Get().theTarget));
289// G4cout << aHadron.GetMomentum()<<" ";
290// G4cout << aHadron.GetTotalMomentum()<<G4endl;
291}
static constexpr double twopi
Definition: G4SIunits.hh:56
#define G4UniformRand()
Definition: Randomize.hh:52
G4double SampleMax(G4double energy)
G4double Sample(G4double x)
void SetMomentum(const G4double x, const G4double y, const G4double z)
void SetTotalEnergy(const G4double en)
G4double GetKineticEnergy() const
G4double GetTotalEnergy() const
G4ThreeVector GetMomentum() const
void Lorentz(const G4ReactionProduct &p1, const G4ReactionProduct &p2)
G4double GetMass() const
void SetMass(const G4double mas)

References fCache, frameFlag, G4cout, G4endl, G4UniformRand, G4ReactionProduct::GetKineticEnergy(), G4ReactionProduct::GetMass(), G4ReactionProduct::GetMomentum(), G4ReactionProduct::GetTotalEnergy(), G4ReactionProduct::GetTotalMomentum(), G4ReactionProduct::Lorentz(), G4ParticleHPPartial::Sample(), G4ParticleHPLegendreStore::SampleMax(), G4ReactionProduct::SetKineticEnergy(), G4ReactionProduct::SetMass(), G4ReactionProduct::SetMomentum(), G4ReactionProduct::SetTotalEnergy(), theAngularDistributionType, theCoefficients, theIsoFlag, theProbArray, and twopi.

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

◆ SetProjectileRP()

void G4ParticleHPAngular::SetProjectileRP ( const G4ReactionProduct anIncidentParticleRP)
inline

◆ SetTarget()

void G4ParticleHPAngular::SetTarget ( const G4ReactionProduct aTarget)
inline

Field Documentation

◆ fCache

G4Cache<toBeCached> G4ParticleHPAngular::fCache
private

◆ frameFlag

G4int G4ParticleHPAngular::frameFlag
private

Definition at line 103 of file G4ParticleHPAngular.hh.

Referenced by G4ParticleHPAngular(), Init(), and SampleAndUpdate().

◆ targetMass

G4double G4ParticleHPAngular::targetMass
private

Definition at line 112 of file G4ParticleHPAngular.hh.

Referenced by G4ParticleHPAngular(), GetTargetMass(), and Init().

◆ theAngularDistributionType

G4int G4ParticleHPAngular::theAngularDistributionType
private

Definition at line 102 of file G4ParticleHPAngular.hh.

Referenced by G4ParticleHPAngular(), Init(), and SampleAndUpdate().

◆ theCoefficients

G4ParticleHPLegendreStore* G4ParticleHPAngular::theCoefficients
private

◆ theIsoFlag

G4bool G4ParticleHPAngular::theIsoFlag
private

Definition at line 105 of file G4ParticleHPAngular.hh.

Referenced by G4ParticleHPAngular(), Init(), and SampleAndUpdate().

◆ theProbArray

G4ParticleHPPartial* G4ParticleHPAngular::theProbArray
private

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