| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // This software was developed by Lawrence Liv 27 // This software was developed by Lawrence Livermore National Laboratory.
28 // 28 //
29 // Redistribution and use in source and binary 29 // Redistribution and use in source and binary forms, with or without
30 // modification, are permitted provided that t 30 // modification, are permitted provided that the following conditions are met:
31 // 31 //
32 // 1. Redistributions of source code must reta 32 // 1. Redistributions of source code must retain the above copyright notice,
33 // this list of conditions and the following 33 // this list of conditions and the following disclaimer.
34 // 2. Redistributions in binary form must repr 34 // 2. Redistributions in binary form must reproduce the above copyright notice,
35 // this list of conditions and the following 35 // this list of conditions and the following disclaimer in the documentation
36 // and/or other materials provided with the 36 // and/or other materials provided with the distribution.
37 // 3. The name of the author may not be used t 37 // 3. The name of the author may not be used to endorse or promote products
38 // derived from this software without specif 38 // derived from this software without specific prior written permission.
39 // 39 //
40 // THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``A 40 // THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR IMPLIED
41 // WARRANTIES, INCLUDING, BUT NOT LIMITED TO, 41 // WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
42 // MERCHANTABILITY AND FITNESS FOR A PARTICULA 42 // MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO
43 // EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DI 43 // EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
44 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGE 44 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
45 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES 45 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
46 // OR BUSINESS INTERRUPTION) HOWEVER CAUSED AN 46 // OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
47 // WHETHER IN CONTRACT, STRICT LIABILITY, OR T 47 // WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
48 // OTHERWISE) ARISING IN ANY WAY OUT OF THE US 48 // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
49 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 49 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
50 // 50 //
51 // Copyright (c) 2006 The Regents of the Unive 51 // Copyright (c) 2006 The Regents of the University of California.
52 // All rights reserved. 52 // All rights reserved.
53 // UCRL-CODE-224807 53 // UCRL-CODE-224807
54 // 54 //
55 // 55 //
56 // 56 //
57 // neutron_hp -- source file 57 // neutron_hp -- source file
58 // J.M. Verbeke, Jan-2007 58 // J.M. Verbeke, Jan-2007
59 // A low energy neutron-induced fission model. 59 // A low energy neutron-induced fission model.
60 // 60 //
61 61
62 #include "G4FissionLibrary.hh" 62 #include "G4FissionLibrary.hh"
63 #include "G4ParticleHPManager.hh" 63 #include "G4ParticleHPManager.hh"
64 #include "G4SystemOfUnits.hh" 64 #include "G4SystemOfUnits.hh"
65 #include "G4PhysicsModelCatalog.hh" 65 #include "G4PhysicsModelCatalog.hh"
66 66
67 G4FissionLibrary::G4FissionLibrary() 67 G4FissionLibrary::G4FissionLibrary()
68 : G4ParticleHPFinalState(), theIsotope(0), t 68 : G4ParticleHPFinalState(), theIsotope(0), targetMass(0.0), secID(-1)
69 { 69 {
70 hasXsec = false; 70 hasXsec = false;
71 fe=0; 71 fe=0;
72 secID = G4PhysicsModelCatalog::GetModelID( " 72 secID = G4PhysicsModelCatalog::GetModelID( "model_G4LLNLFission" );
73 } 73 }
74 74
75 G4FissionLibrary::~G4FissionLibrary() 75 G4FissionLibrary::~G4FissionLibrary()
76 {} 76 {}
77 77
78 G4ParticleHPFinalState * G4FissionLibrary::New 78 G4ParticleHPFinalState * G4FissionLibrary::New()
79 { 79 {
80 G4FissionLibrary * theNew = new G4FissionLib 80 G4FissionLibrary * theNew = new G4FissionLibrary;
81 return theNew; 81 return theNew;
82 } 82 }
83 83
84 //void G4FissionLibrary::Init (G4double A, G4d 84 //void G4FissionLibrary::Init (G4double A, G4double Z, G4String & dirName, G4String &)
85 void G4FissionLibrary::Init (G4double A, G4dou << 85 void G4FissionLibrary::Init (G4double A, G4double Z, G4int M, G4String & dirName, G4String &, G4ParticleDefinition*)
86 { 86 {
87 G4String tString = "/FS/"; 87 G4String tString = "/FS/";
88 G4bool dbool; 88 G4bool dbool;
89 theIsotope = static_cast<G4int>(1000*Z+A); 89 theIsotope = static_cast<G4int>(1000*Z+A);
90 G4ParticleHPDataUsed aFile = theNames.GetNam 90 G4ParticleHPDataUsed aFile = theNames.GetName(static_cast<G4int>(A), static_cast<G4int>(Z), M, dirName, tString, dbool);
91 G4String filename = aFile.GetName(); 91 G4String filename = aFile.GetName();
92 92
93 if(!dbool) 93 if(!dbool)
94 { 94 {
95 hasAnyData = false; 95 hasAnyData = false;
96 hasFSData = false; 96 hasFSData = false;
97 hasXsec = false; 97 hasXsec = false;
98 return; 98 return;
99 } 99 }
100 //std::ifstream theData(filename, std::ios:: 100 //std::ifstream theData(filename, std::ios::in);
101 std::istringstream theData(std::ios::in); 101 std::istringstream theData(std::ios::in);
102 G4ParticleHPManager::GetInstance()->GetDataS 102 G4ParticleHPManager::GetInstance()->GetDataStream(filename,theData);
103 103
104 // here it comes 104 // here it comes
105 G4int infoType, dataType; 105 G4int infoType, dataType;
106 hasFSData = false; 106 hasFSData = false;
107 while (theData >> infoType) // Loop checking 107 while (theData >> infoType) // Loop checking, 11.03.2015, T. Koi
108 { 108 {
109 hasFSData = true; 109 hasFSData = true;
110 theData >> dataType; 110 theData >> dataType;
111 switch(infoType) 111 switch(infoType)
112 { 112 {
113 case 1: 113 case 1:
114 if(dataType==4) theNeutronAngularDis.I 114 if(dataType==4) theNeutronAngularDis.Init(theData);
115 if(dataType==5) thePromptNeutronEnDis. 115 if(dataType==5) thePromptNeutronEnDis.Init(theData);
116 if(dataType==12) theFinalStatePhotons. 116 if(dataType==12) theFinalStatePhotons.InitMean(theData);
117 if(dataType==14) theFinalStatePhotons. 117 if(dataType==14) theFinalStatePhotons.InitAngular(theData);
118 if(dataType==15) theFinalStatePhotons. 118 if(dataType==15) theFinalStatePhotons.InitEnergies(theData);
119 break; 119 break;
120 case 2: 120 case 2:
121 if(dataType==1) theFinalStateNeutrons. 121 if(dataType==1) theFinalStateNeutrons.InitMean(theData);
122 break; 122 break;
123 case 3: 123 case 3:
124 if(dataType==1) theFinalStateNeutrons. 124 if(dataType==1) theFinalStateNeutrons.InitDelayed(theData);
125 if(dataType==5) theDelayedNeutronEnDis 125 if(dataType==5) theDelayedNeutronEnDis.Init(theData);
126 break; 126 break;
127 case 4: 127 case 4:
128 if(dataType==1) theFinalStateNeutrons. 128 if(dataType==1) theFinalStateNeutrons.InitPrompt(theData);
129 break; 129 break;
130 case 5: 130 case 5:
131 if(dataType==1) theEnergyRelease.Init( 131 if(dataType==1) theEnergyRelease.Init(theData);
132 break; 132 break;
133 default: 133 default:
134 G4cout << "G4FissionLibrary::Init: unk 134 G4cout << "G4FissionLibrary::Init: unknown data type"<<dataType<<G4endl;
135 throw G4HadronicException(__FILE__, __ 135 throw G4HadronicException(__FILE__, __LINE__, "G4FissionLibrary::Init: unknown data type");
136 break; 136 break;
137 } 137 }
138 } 138 }
139 targetMass = theFinalStateNeutrons.GetTarget 139 targetMass = theFinalStateNeutrons.GetTargetMass();
140 //theData.close(); 140 //theData.close();
141 } 141 }
142 142
143 G4HadFinalState* G4FissionLibrary::ApplyYourse 143 G4HadFinalState* G4FissionLibrary::ApplyYourself(const G4HadProjectile & theTrack)
144 { 144 {
145 145
146 if ( theResult.Get() == NULL ) theResult.Put 146 if ( theResult.Get() == NULL ) theResult.Put( new G4HadFinalState );
147 theResult.Get()->Clear(); 147 theResult.Get()->Clear();
148 148
149 // prepare neutron 149 // prepare neutron
150 G4double eKinetic = theTrack.GetKineticEnerg 150 G4double eKinetic = theTrack.GetKineticEnergy();
151 const G4HadProjectile* incidentParticle = &t 151 const G4HadProjectile* incidentParticle = &theTrack;
152 G4ReactionProduct theNeutron(incidentParticl 152 G4ReactionProduct theNeutron(incidentParticle->GetDefinition() );
153 theNeutron.SetMomentum(incidentParticle->Get 153 theNeutron.SetMomentum(incidentParticle->Get4Momentum().vect() );
154 theNeutron.SetKineticEnergy(eKinetic); 154 theNeutron.SetKineticEnergy(eKinetic);
155 155
156 // prepare target 156 // prepare target
157 G4Nucleus aNucleus; 157 G4Nucleus aNucleus;
158 G4ReactionProduct theTarget; 158 G4ReactionProduct theTarget;
159 G4ThreeVector neuVelo = (1./incidentParticle 159 G4ThreeVector neuVelo = (1./incidentParticle->GetDefinition()->GetPDGMass())*theNeutron.GetMomentum();
160 theTarget = aNucleus.GetBiasedThermalNucleus 160 theTarget = aNucleus.GetBiasedThermalNucleus( targetMass, neuVelo, theTrack.GetMaterial()->GetTemperature());
161 161
162 // set neutron and target in the FS classes 162 // set neutron and target in the FS classes
163 //theNeutronAngularDis.SetNeutron(theNeutron 163 //theNeutronAngularDis.SetNeutron(theNeutron);
164 theNeutronAngularDis.SetProjectileRP(theNeut 164 theNeutronAngularDis.SetProjectileRP(theNeutron);
165 theNeutronAngularDis.SetTarget(theTarget); 165 theNeutronAngularDis.SetTarget(theTarget);
166 166
167 // boost to target rest system 167 // boost to target rest system
168 theNeutron.Lorentz(theNeutron, -1*theTarget) 168 theNeutron.Lorentz(theNeutron, -1*theTarget);
169 169
170 eKinetic = theNeutron.GetKineticEnergy(); 170 eKinetic = theNeutron.GetKineticEnergy();
171 171
172 // dice neutron and gamma multiplicities, en 172 // dice neutron and gamma multiplicities, energies and momenta in Lab. @@
173 // no energy conservation on an event-to-eve 173 // no energy conservation on an event-to-event basis. we rely on the data to be ok. @@
174 // also for mean, we rely on the consistency 174 // also for mean, we rely on the consistency of the data. @@
175 175
176 G4int nPrompt=0, gPrompt=0; 176 G4int nPrompt=0, gPrompt=0;
177 SampleMult(theTrack, &nPrompt, &gPrompt, eKi 177 SampleMult(theTrack, &nPrompt, &gPrompt, eKinetic);
178 178
179 // Build neutrons and add them to dynamic pa 179 // Build neutrons and add them to dynamic particle vector
180 G4double momentum; 180 G4double momentum;
181 for(G4int i=0; i<nPrompt; i++) 181 for(G4int i=0; i<nPrompt; i++)
182 { 182 {
183 G4DynamicParticle * it = new G4DynamicPart 183 G4DynamicParticle * it = new G4DynamicParticle;
184 it->SetDefinition(G4Neutron::Neutron()); 184 it->SetDefinition(G4Neutron::Neutron());
185 it->SetKineticEnergy(fe->getNeutronEnergy( 185 it->SetKineticEnergy(fe->getNeutronEnergy(i)*MeV);
186 momentum = it->GetTotalMomentum(); 186 momentum = it->GetTotalMomentum();
187 G4ThreeVector temp(momentum*fe->getNeutron 187 G4ThreeVector temp(momentum*fe->getNeutronDircosu(i),
188 momentum*fe->getNeutron 188 momentum*fe->getNeutronDircosv(i),
189 momentum*fe->getNeutron 189 momentum*fe->getNeutronDircosw(i));
190 it->SetMomentum( temp ); 190 it->SetMomentum( temp );
191 // it->SetGlobalTime(fe->getNeutronAge(i)*s 191 // it->SetGlobalTime(fe->getNeutronAge(i)*second);
192 theResult.Get()->AddSecondary(it, secID); 192 theResult.Get()->AddSecondary(it, secID);
193 // G4cout <<"G4FissionLibrary::ApplyYoursel 193 // G4cout <<"G4FissionLibrary::ApplyYourself: energy of prompt neutron " << i << " = " << it->GetKineticEnergy()<<G4endl;
194 } 194 }
195 195
196 // Build gammas, lorentz transform them, and 196 // Build gammas, lorentz transform them, and add them to dynamic particle vector
197 for(G4int i=0; i<gPrompt; i++) 197 for(G4int i=0; i<gPrompt; i++)
198 { 198 {
199 G4ReactionProduct * thePhoton = new G4Reac 199 G4ReactionProduct * thePhoton = new G4ReactionProduct;
200 thePhoton->SetDefinition(G4Gamma::Gamma()) 200 thePhoton->SetDefinition(G4Gamma::Gamma());
201 thePhoton->SetKineticEnergy(fe->getPhotonE 201 thePhoton->SetKineticEnergy(fe->getPhotonEnergy(i)*MeV);
202 momentum = thePhoton->GetTotalMomentum(); 202 momentum = thePhoton->GetTotalMomentum();
203 G4ThreeVector temp(momentum*fe->getPhotonD 203 G4ThreeVector temp(momentum*fe->getPhotonDircosu(i),
204 momentum*fe->getPhotonD 204 momentum*fe->getPhotonDircosv(i),
205 momentum*fe->getPhotonD 205 momentum*fe->getPhotonDircosw(i));
206 thePhoton->SetMomentum( temp ); 206 thePhoton->SetMomentum( temp );
207 thePhoton->Lorentz(*thePhoton, -1.*theTarg 207 thePhoton->Lorentz(*thePhoton, -1.*theTarget);
208 208
209 G4DynamicParticle * it = new G4DynamicPart 209 G4DynamicParticle * it = new G4DynamicParticle;
210 it->SetDefinition(thePhoton->GetDefinition 210 it->SetDefinition(thePhoton->GetDefinition());
211 it->SetMomentum(thePhoton->GetMomentum()); 211 it->SetMomentum(thePhoton->GetMomentum());
212 // it->SetGlobalTime(fe->getPhotonAge(i)*se 212 // it->SetGlobalTime(fe->getPhotonAge(i)*second);
213 // G4cout <<"G4FissionLibrary::ApplyYoursel 213 // G4cout <<"G4FissionLibrary::ApplyYourself: energy of prompt photon " << i << " = " << it->GetKineticEnergy()<<G4endl;
214 theResult.Get()->AddSecondary(it, secID); 214 theResult.Get()->AddSecondary(it, secID);
215 delete thePhoton; 215 delete thePhoton;
216 } 216 }
217 // G4cout <<"G4FissionLibrary::ApplyYourself: 217 // G4cout <<"G4FissionLibrary::ApplyYourself: Number of secondaries = "<<theResult.GetNumberOfSecondaries()<< G4endl;
218 // G4cout <<"G4FissionLibrary::ApplyYourself: 218 // G4cout <<"G4FissionLibrary::ApplyYourself: Number of induced prompt neutron = "<<nPrompt<<G4endl;
219 // G4cout <<"G4FissionLibrary::ApplyYourself: 219 // G4cout <<"G4FissionLibrary::ApplyYourself: Number of induced prompt photons = "<<gPrompt<<G4endl;
220 220
221 // finally deal with local energy deposition 221 // finally deal with local energy depositions.
222 G4double eDepByFragments = theEnergyRelease. 222 G4double eDepByFragments = theEnergyRelease.GetFragmentKinetic();
223 theResult.Get()->SetLocalEnergyDeposit(eDepB 223 theResult.Get()->SetLocalEnergyDeposit(eDepByFragments);
224 // G4cout << "G4FissionLibrary::local energy 224 // G4cout << "G4FissionLibrary::local energy deposit" << eDepByFragments<<G4endl;
225 // clean up the primary neutron 225 // clean up the primary neutron
226 theResult.Get()->SetStatusChange(stopAndKill 226 theResult.Get()->SetStatusChange(stopAndKill);
227 return theResult.Get(); 227 return theResult.Get();
228 } 228 }
229 229
230 void G4FissionLibrary::SampleMult(const G4HadP 230 void G4FissionLibrary::SampleMult(const G4HadProjectile & theTrack, G4int* nPrompt,
231 G4int* gPro 231 G4int* gPrompt, G4double eKinetic)
232 { 232 {
233 G4double promptNeutronMulti = 0; 233 G4double promptNeutronMulti = 0;
234 promptNeutronMulti = theFinalStateNeutrons. 234 promptNeutronMulti = theFinalStateNeutrons.GetPrompt(eKinetic); // prompt nubar from Geant
235 G4double delayedNeutronMulti = 0; 235 G4double delayedNeutronMulti = 0;
236 delayedNeutronMulti = theFinalStateNeutrons 236 delayedNeutronMulti = theFinalStateNeutrons.GetDelayed(eKinetic); // delayed nubar from Geant
237 237
238 G4double time = theTrack.GetGlobalTime()/se 238 G4double time = theTrack.GetGlobalTime()/second;
239 G4double totalNeutronMulti = theFinalStateN 239 G4double totalNeutronMulti = theFinalStateNeutrons.GetMean(eKinetic);
240 if(delayedNeutronMulti==0&&promptNeutronMul 240 if(delayedNeutronMulti==0&&promptNeutronMulti==0) {
241 // no data for prompt and delayed neutron 241 // no data for prompt and delayed neutrons in Geant
242 // but there is perhaps data for the tota 242 // but there is perhaps data for the total neutron multiplicity, in which case
243 // we use it for prompt neutron emission 243 // we use it for prompt neutron emission
244 if (fe != 0) delete fe; 244 if (fe != 0) delete fe;
245 fe = new G4fissionEvent(theIsotope, time, 245 fe = new G4fissionEvent(theIsotope, time, totalNeutronMulti, eKinetic);
246 } else { 246 } else {
247 // prompt nubar != 0 || delayed nubar != 247 // prompt nubar != 0 || delayed nubar != 0
248 if (fe != 0) delete fe; 248 if (fe != 0) delete fe;
249 fe = new G4fissionEvent(theIsotope, time, 249 fe = new G4fissionEvent(theIsotope, time, promptNeutronMulti, eKinetic);
250 } 250 }
251 *nPrompt = fe->getNeutronNu(); 251 *nPrompt = fe->getNeutronNu();
252 if (*nPrompt == -1) *nPrompt = 0; // the fi 252 if (*nPrompt == -1) *nPrompt = 0; // the fission library libFission.a has no data for neutrons
253 *gPrompt = fe->getPhotonNu(); 253 *gPrompt = fe->getPhotonNu();
254 if (*gPrompt == -1) *gPrompt = 0; // the fi 254 if (*gPrompt == -1) *gPrompt = 0; // the fission library libFission.a has no data for gammas
255 } 255 }
256 256
257 257