| 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 // Class Description 26 // Class Description
27 // Manager of NetronHP 27 // Manager of NetronHP
28 // << 28 //
29 // 121031 First implementation done by T. Koi 29 // 121031 First implementation done by T. Koi (SLAC/PPA)
30 // P. Arce, June-2014 Conversion neutron_hp to 30 // P. Arce, June-2014 Conversion neutron_hp to particle_hp
31 // V. Ivanchenko, July-2023 Basic revision of <<
32 // 31 //
33 #include "G4ParticleHPManager.hh" 32 #include "G4ParticleHPManager.hh"
34 <<
35 #include "G4Exception.hh" <<
36 #include "G4HadronicException.hh" <<
37 #include "G4ParticleHPMessenger.hh" <<
38 #include "G4ParticleDefinition.hh" <<
39 #include "G4HadronicParameters.hh" <<
40 #include "G4ParticleHPThreadLocalManager.hh" 33 #include "G4ParticleHPThreadLocalManager.hh"
41 #include "G4SystemOfUnits.hh" << 34 #include "G4ParticleHPMessenger.hh"
42 << 35 #include "G4HadronicException.hh"
43 #include <zlib.h> <<
44 #include <fstream> <<
45 36
46 G4ParticleHPManager* G4ParticleHPManager::inst << 37 //G4ThreadLocal G4ParticleHPManager* G4ParticleHPManager::instance = NULL;
>> 38 G4ParticleHPManager* G4ParticleHPManager::instance = G4ParticleHPManager::GetInstance();
47 39
48 G4ParticleHPManager::G4ParticleHPManager() 40 G4ParticleHPManager::G4ParticleHPManager()
49 : theMinEnergyDBRC(0.1 * CLHEP::eV), << 41 :/*RWB(NULL),*/
50 theMaxEnergyDBRC(210. * CLHEP::eV), << 42 verboseLevel(1)
51 theMaxEnergyDoppler(30. * CLHEP::keV) << 43 ,USE_ONLY_PHOTONEVAPORATION(false)
52 { << 44 ,SKIP_MISSING_ISOTOPES(false)
53 messenger = new G4ParticleHPMessenger(this); << 45 ,NEGLECT_DOPPLER(false)
54 verboseLevel = G4HadronicParameters::Instanc << 46 ,DO_NOT_ADJUST_FINAL_STATE(false)
55 const char* ss = G4FindDataDir("NeutronHPNam << 47 ,PRODUCE_FISSION_FRAGMENTS(false)
56 if (nullptr != ss) { CHECK_HP_NAMES = true; << 48 ,USE_NRESP71_MODEL(false)
57 ss = G4FindDataDir("G4PHP_DO_NOT_CHECK_DIFF_ << 49 ,theElasticCrossSections(NULL)
58 if (nullptr != ss) { PHP_CHECK = false; } << 50 ,theCaptureCrossSections(NULL)
59 ss = G4FindDataDir("G4PHP_MULTIPLICITY_METHO << 51 //,theInelasticCrossSections(NULL)
60 if (nullptr != ss && "BetweenInts" == G4Stri << 52 ,theFissionCrossSections(NULL)
61 ss = G4FindDataDir("G4ParticleHPDebug"); << 53 ,theElasticFSs(NULL)
62 if (nullptr != ss) { DEBUG = true; } << 54 //,theInelasticFSs(NULL)
63 << 55 ,theCaptureFSs(NULL)
64 // identify and check data path once - it sh << 56 ,theFissionFSs(NULL)
65 const char* nch = G4FindDataDir("G4NEUTRONHP << 57 ,theTSCoherentCrossSections(NULL)
66 if (nullptr == nch) { << 58 ,theTSIncoherentCrossSections(NULL)
67 G4Exception("G4ParticleHPManager::G4Partic << 59 ,theTSInelasticCrossSections(NULL)
68 FatalException, "G4NEUTRONXSDA << 60 ,theTSCoherentFinalStates(NULL)
69 } else { << 61 ,theTSIncoherentFinalStates(NULL)
70 fDataPath[0] = G4String(nch); << 62 ,theTSInelasticFinalStates(NULL)
71 } << 63 {
72 // path may be defined by two environment va << 64 messenger = new G4ParticleHPMessenger( this );
73 // it is not mandatory to access PHP data - << 65 if ( getenv( "G4NEUTRONHP_DO_NOT_ADJUST_FINAL_STATE" ) || getenv("G4PHP_DO_NOT_ADJUST_FINAL_STATE") ) DO_NOT_ADJUST_FINAL_STATE = true;
74 const char* ttp = G4FindDataDir("G4PARTICLEH << 66 if ( getenv( "G4NEUTRONHP_USE_ONLY_PHOTONEVAPORATION" ) ) USE_ONLY_PHOTONEVAPORATION = true;
75 G4String tendl = (nullptr == ttp) ? G4String << 67 if ( getenv( "G4NEUTRONHP_NEGLECT_DOPPLER" ) || getenv("G4PHP_NEGLECT_DOPPLER") ) NEGLECT_DOPPLER = true;
76 const char* ssp = G4FindDataDir("G4PROTONHPD << 68 if ( getenv( "G4NEUTRONHP_SKIP_MISSING_ISOTOPES" ) ) SKIP_MISSING_ISOTOPES = true;
77 fDataPath[1] = (nullptr == ssp) ? tendl + "/ << 69 if ( getenv( "G4NEUTRONHP_PRODUCE_FISSION_FRAGMENTS" ) ) PRODUCE_FISSION_FRAGMENTS = true;
78 << 70 if ( getenv( "G4PHP_USE_NRESP71_MODEL" ) ) USE_NRESP71_MODEL = true;
79 ssp = G4FindDataDir("G4DEUTERONHPDATA"); <<
80 fDataPath[2] = (nullptr == ssp) ? tendl + "/ <<
81 <<
82 ssp = G4FindDataDir("G4TRITONHPDATA"); <<
83 fDataPath[3] = (nullptr == ssp) ? tendl + "/ <<
84 <<
85 ssp = G4FindDataDir("G4HE3HPDATA"); <<
86 fDataPath[4] = (nullptr == ssp) ? tendl + "/ <<
87 <<
88 ssp = G4FindDataDir("G4ALPHAHPDATA"); <<
89 fDataPath[5] = (nullptr == ssp) ? tendl + "/ <<
90 } 71 }
91 <<
92 G4ParticleHPManager::~G4ParticleHPManager() 72 G4ParticleHPManager::~G4ParticleHPManager()
93 { 73 {
94 delete messenger; << 74 delete messenger;
95 } 75 }
96 <<
97 G4ParticleHPManager* G4ParticleHPManager::GetI <<
98 { <<
99 static G4ParticleHPManager manager; <<
100 if (instance == nullptr) { <<
101 instance = &manager; <<
102 } <<
103 return instance; <<
104 } <<
105 <<
106 G4int G4ParticleHPManager::GetPHPIndex(const G <<
107 G4int pdg = part->GetPDGEncoding(); <<
108 G4int idx; <<
109 if (pdg == 2112) { idx = 0; } <<
110 else if (pdg == 2212) { idx = 1; } <<
111 else if (pdg == 1000010020) { idx = 2; } <<
112 else if (pdg == 1000010030) { idx = 3; } <<
113 else if (pdg == 1000020030) { idx = 4; } <<
114 else if (pdg == 1000020040) { idx = 5; } <<
115 else { <<
116 idx = 0; <<
117 G4ExceptionDescription ed; <<
118 ed << "Particle " << part->GetParticleName <<
119 << " cannot be handled by the ParticleH <<
120 G4Exception("G4ParticleHPManager::G4Partic <<
121 FatalException, ed, ""); <<
122 } <<
123 return idx; <<
124 } <<
125 <<
126 const G4String& <<
127 G4ParticleHPManager::GetParticleHPPath(const G <<
128 return fDataPath[GetPHPIndex(part)]; <<
129 } <<
130 <<
131 void G4ParticleHPManager::OpenReactionWhiteBoa 76 void G4ParticleHPManager::OpenReactionWhiteBoard()
132 { 77 {
133 G4ParticleHPThreadLocalManager::GetInstance( << 78 // if ( RWB != NULL ) {
>> 79 // G4cout << "Warning: G4ParticleHPReactionWhiteBoard is tried doubly opening" << G4endl;
>> 80 // RWB = new G4ParticleHPReactionWhiteBoard();
>> 81 // }
>> 82 //
>> 83 // RWB = new G4ParticleHPReactionWhiteBoard();
>> 84 G4ParticleHPThreadLocalManager::GetInstance()->OpenReactionWhiteBoard();
134 } 85 }
135 <<
136 G4ParticleHPReactionWhiteBoard* G4ParticleHPMa 86 G4ParticleHPReactionWhiteBoard* G4ParticleHPManager::GetReactionWhiteBoard()
137 { 87 {
138 return G4ParticleHPThreadLocalManager::GetIn << 88 // if ( RWB == NULL ) {
>> 89 // G4cout << "Warning: try to access G4ParticleHPReactionWhiteBoard before opening" << G4endl;
>> 90 // RWB = new G4ParticleHPReactionWhiteBoard();
>> 91 // }
>> 92 // return RWB;
>> 93 return G4ParticleHPThreadLocalManager::GetInstance()->GetReactionWhiteBoard();
139 } 94 }
140 <<
141 void G4ParticleHPManager::CloseReactionWhiteBo 95 void G4ParticleHPManager::CloseReactionWhiteBoard()
142 { 96 {
143 G4ParticleHPThreadLocalManager::GetInstance( << 97 G4ParticleHPThreadLocalManager::GetInstance()->CloseReactionWhiteBoard();
144 } 98 }
145 99
146 void G4ParticleHPManager::GetDataStream(const << 100 #include "zlib.h"
>> 101 #include <fstream>
>> 102 void G4ParticleHPManager::GetDataStream( G4String filename , std::istringstream& iss )
147 { 103 {
148 G4String* data = nullptr; << 104 //if ( getenv( "TEST04" ) && filename != "INVALID" ) G4cout << "Reading " << filename << G4endl;
149 G4String compfilename(filename); << 105 G4String* data=NULL;
150 compfilename += ".z"; << 106 G4String compfilename(filename);
151 auto in = new std::ifstream(compfilename, st << 107 compfilename += ".z";
152 if (in->good()) { << 108 std::ifstream* in = new std::ifstream ( compfilename , std::ios::binary | std::ios::ate );
153 // Use the compressed file << 109 if ( in->good() )
154 std::streamoff file_size = in->tellg(); << 110 {
155 in->seekg(0, std::ios::beg); << 111 // Use the compressed file
156 auto compdata = new Bytef[file_size]; << 112 G4int file_size = in->tellg();
157 << 113 in->seekg( 0 , std::ios::beg );
158 while (*in) { // Loop checking, 11.05.201 << 114 Bytef* compdata = new Bytef[ file_size ];
159 in->read((char*)compdata, file_size); << 115
160 } << 116 while ( *in ) { // Loop checking, 11.05.2015, T. Koi
161 << 117 in->read( (char*)compdata , file_size );
162 auto complen = (uLongf)(file_size * 4); <<
163 auto uncompdata = new Bytef[complen]; <<
164 <<
165 while (Z_OK != uncompress(uncompdata, &com <<
166 { // Loop checking, 11.05.2015, T. Koi <<
167 delete[] uncompdata; <<
168 complen *= 2; <<
169 uncompdata = new Bytef[complen]; <<
170 } <<
171 delete[] compdata; <<
172 // Now "complen" has uncomplessed size <<
173 data = new G4String((char*)uncompdata, (G4 <<
174 delete[] uncompdata; <<
175 } <<
176 else { <<
177 // Use regular text file <<
178 std::ifstream thefData(filename, std::ios: <<
179 if (thefData.good()) { <<
180 std::streamoff file_size = thefData.tell <<
181 thefData.seekg(0, std::ios::beg); <<
182 auto filedata = new char[file_size]; <<
183 while (thefData) { // Loop checking, 11 <<
184 thefData.read(filedata, file_size); <<
185 } 118 }
186 thefData.close(); <<
187 data = new G4String(filedata, file_size) <<
188 delete[] filedata; <<
189 } <<
190 else { <<
191 // found no data file <<
192 // set error bit to the stream <<
193 iss.setstate(std::ios::badbit); <<
194 } <<
195 } <<
196 if (data != nullptr) { <<
197 iss.str(*data); <<
198 G4String id; <<
199 iss >> id; <<
200 if (id == "G4NDL") { <<
201 // Register information of file <<
202 G4String source; <<
203 iss >> source; <<
204 register_data_file(filename, source); <<
205 } <<
206 else { <<
207 iss.seekg(0, std::ios::beg); <<
208 } <<
209 } <<
210 in->close(); <<
211 delete in; <<
212 delete data; <<
213 } <<
214 119
215 void G4ParticleHPManager::GetDataStream2(const << 120 uLongf complen = (uLongf) ( file_size*4 );
216 { << 121 Bytef* uncompdata = new Bytef[complen];
217 // Checking existance of data file <<
218 122
219 G4String compfilename(filename); << 123 while ( Z_OK != uncompress ( uncompdata , &complen , compdata , file_size ) ) { // Loop checking, 11.05.2015, T. Koi
220 compfilename += ".z"; << 124 //G4cout << "Too small, retry 2 times bigger size." << G4endl;
221 auto in = new std::ifstream(compfilename, st << 125 delete[] uncompdata;
222 if (in->good()) { << 126 complen *= 2;
223 // Compressed file is exist << 127 uncompdata = new Bytef[complen];
224 in->close(); << 128 }
225 } << 129 delete [] compdata;
226 else { << 130 // Now "complen" has uncomplessed size
227 std::ifstream thefData(filename, std::ios: << 131 data = new G4String ( (char*)uncompdata , (G4long)complen );
228 if (thefData.good()) { << 132 delete [] uncompdata;
229 // Regular text file is exist << 133 } else {
230 thefData.close(); << 134 // Use regular text file
231 } << 135 std::ifstream thefData( filename , std::ios::in | std::ios::ate );
232 else { << 136 if ( thefData.good() ) {
233 // found no data file << 137 G4int file_size = thefData.tellg();
234 // set error bit to the stream << 138 thefData.seekg( 0 , std::ios::beg );
235 iss.setstate(std::ios::badbit); << 139 char* filedata = new char[ file_size ];
236 } << 140 while ( thefData ) { // Loop checking, 11.05.2015, T. Koi
237 } << 141 thefData.read( filedata , file_size );
238 delete in; << 142 }
>> 143 thefData.close();
>> 144 data = new G4String ( filedata , file_size );
>> 145 delete [] filedata;
>> 146 } else {
>> 147 // found no data file
>> 148 // set error bit to the stream
>> 149 iss.setstate( std::ios::badbit );
>> 150 }
>> 151 }
>> 152 if ( data != NULL ) {
>> 153 iss.str(*data);
>> 154 G4String id;
>> 155 iss >> id;
>> 156 if ( id == "G4NDL" ) {
>> 157 //Register information of file
>> 158 G4String source;
>> 159 iss >> source;
>> 160 register_data_file(filename,source);
>> 161 } else {
>> 162 iss.seekg( 0 , std::ios::beg );
>> 163 }
>> 164 }
>> 165 //G4cout << iss.rdbuf()->in_avail() << G4endl;
>> 166 in->close(); delete in;
>> 167 delete data;
>> 168 }
>> 169 // Checking existance of data file
>> 170 void G4ParticleHPManager::GetDataStream2( G4String filename , std::istringstream& iss )
>> 171 {
>> 172 G4String compfilename(filename);
>> 173 compfilename += ".z";
>> 174 std::ifstream* in = new std::ifstream ( compfilename , std::ios::binary | std::ios::ate );
>> 175 if ( in->good() )
>> 176 {
>> 177 // Compressed file is exist
>> 178 in->close();
>> 179 } else {
>> 180 std::ifstream thefData( filename , std::ios::in | std::ios::ate );
>> 181 if ( thefData.good() ) {
>> 182 // Regular text file is exist
>> 183 thefData.close();
>> 184 } else {
>> 185 // found no data file
>> 186 // set error bit to the stream
>> 187 iss.setstate( std::ios::badbit );
>> 188 }
>> 189 }
>> 190 delete in;
239 } 191 }
240 192
241 void G4ParticleHPManager::SetVerboseLevel(G4in << 193 void G4ParticleHPManager::SetVerboseLevel( G4int newValue )
242 { 194 {
243 G4cout << "You are setting a new verbose lev << 195 G4cout << "You are setting a new verbose level for Particle HP package." << G4endl;
244 G4cout << "the new value will be used in who << 196 G4cout << "the new value will be used in whole of the Particle HP package, i.e., models and cross sections for Capture, Elastic, Fission and Inelastic interaction." << G4endl;
245 "cross sections for Capture, Elast << 197 verboseLevel = newValue;
246 << G4endl; <<
247 verboseLevel = newValue; <<
248 } 198 }
249 199
250 void G4ParticleHPManager::register_data_file(c << 200 void G4ParticleHPManager::register_data_file(G4String filename, G4String source)
251 { 201 {
252 mDataEvaluation.insert(std::pair<G4String, G << 202 mDataEvaluation.insert( std::pair < G4String , G4String > ( filename , source ) );
253 } 203 }
254 204
255 void G4ParticleHPManager::DumpDataSource() con << 205 void G4ParticleHPManager::DumpDataSource()
256 { 206 {
257 G4cout << "Data source of this Partile HP ca << 207
258 for (const auto& it : mDataEvaluation) { << 208 G4cout << "Data source of this Partile HP calculation are " << G4endl;
259 G4cout << it.first << " " << it.second << << 209 for ( std::map< G4String , G4String >::iterator
260 } << 210 it = mDataEvaluation.begin() ; it != mDataEvaluation.end() ; it++ ) {
261 G4cout << G4endl; << 211 G4cout << it->first << " " << it->second << G4endl;
>> 212 }
>> 213 G4cout << G4endl;
262 } 214 }
263 215
264 void G4ParticleHPManager::DumpSetting() << 216 G4PhysicsTable* G4ParticleHPManager::GetInelasticCrossSections(const G4ParticleDefinition* particle ){
265 { << 217 if ( theInelasticCrossSections.end() != theInelasticCrossSections.find( particle ) )
266 if(isPrinted) { return; } << 218 return theInelasticCrossSections.find( particle )->second;
267 G4cout << G4endl << 219 else
268 << "================================= << 220 return NULL;
269 << "====== ParticleHP Physics P << 221 }
270 << "================================= << 222
271 << " Use only photo-evaporation << 223 void G4ParticleHPManager::RegisterInelasticCrossSections( const G4ParticleDefinition* particle, G4PhysicsTable* val ){
272 << " Skip missing isotopes << 224 theInelasticCrossSections.insert( std::pair<const G4ParticleDefinition* , G4PhysicsTable* >( particle , val ) );
273 << " Neglect Doppler << 225 }
274 << " Do not adjust final state << 226
275 << " Produce fission fragments << 227 std::vector<G4ParticleHPChannelList*>* G4ParticleHPManager::GetInelasticFinalStates(const G4ParticleDefinition* particle) {
276 << " Use WendtFissionModel << 228 if ( theInelasticFSs.end() != theInelasticFSs.find( particle ) )
277 << " Use NRESP71Model << 229 return theInelasticFSs.find( particle )->second;
278 << " Use DBRC << 230 else
279 << " PHP use Poisson << 231 return NULL;
280 << " PHP check << 232 }
281 << " CHECK HP NAMES << 233
282 << " Enable DEBUG << 234 void G4ParticleHPManager::RegisterInelasticFinalStates( const G4ParticleDefinition* particle , std::vector<G4ParticleHPChannelList*>* val ) {
283 << " Use probability tables from << 235 theInelasticFSs.insert ( std::pair<const G4ParticleDefinition*,std::vector<G4ParticleHPChannelList*>*>( particle , val ) );
284 << "================================= <<
285 isPrinted = true; <<
286 } 236 }
287 237