Geant4  10.03
G4IonStoppingData.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4IonStoppingData.cc 96794 2016-05-09 10:09:30Z gcosmo $
27 //
28 // ===========================================================================
29 // GEANT4 class source file
30 //
31 // Class: G4IonStoppingData
32 //
33 // Base class: G4VIonDEDXTable
34 //
35 // Author: Anton Lechner (Anton.Lechner@cern.ch)
36 //
37 // First implementation: 03. 11. 2009
38 //
39 // Modifications:
40 // 25.10.2010 V.Ivanchenko fixed warnings reported by the Coverity tool
41 // 25.10.2011: new scheme for G4Exception (mma)
42 //
43 //
44 // Class description: Class which can read ion stopping power data from
45 // $G4LEDATA/ion_stopping_data
46 //
47 // Comments:
48 //
49 // ===========================================================================
50 //
51 
52 #include <fstream>
53 #include <sstream>
54 #include <iomanip>
55 
56 #include "G4IonStoppingData.hh"
57 #include "G4PhysicsVector.hh"
58 #include "G4LPhysicsFreeVector.hh"
59 #include "G4PhysicalConstants.hh"
60 #include "G4SystemOfUnits.hh"
61 
62 // #########################################################################
63 
65  subDir( leDirectory ) {
66 
67 }
68 
69 // #########################################################################
70 
72 
73  ClearTable();
74 }
75 
76 // #########################################################################
77 
79  G4int atomicNumberIon, // Atomic number of ion
80  G4int atomicNumberElem // Atomic number of elemental material
81  )
82 {
83  G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
84 
85  G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
86 
87  return (iter == dedxMapElements.end()) ? false : true;
88 }
89 
90 // #########################################################################
91 
93  G4int atomicNumberIon, // Atomic number of ion
94  const G4String& matIdentifier // Name or chemical formula of material
95  )
96 {
97  G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
98 
99  G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
100 
101  return (iter == dedxMapMaterials.end()) ? false : true;
102 }
103 
104 // #########################################################################
105 
107  G4int atomicNumberIon, // Atomic number of ion
108  G4int atomicNumberElem // Atomic number of elemental material
109  )
110 {
111  G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
112 
113  G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
114 
115  return (iter != dedxMapElements.end()) ? iter->second : nullptr;
116 }
117 
118 // #########################################################################
119 
121  G4int atomicNumberIon, // Atomic number of ion
122  const G4String& matIdentifier // Name or chemical formula of material
123  )
124 {
125  G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
126 
127  G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
128 
129  return (iter != dedxMapMaterials.end()) ? iter->second : nullptr;
130 }
131 
132 // #########################################################################
133 
135  G4double kinEnergyPerNucleon, // Kinetic energy per nucleon
136  G4int atomicNumberIon, // Atomic number of ion
137  G4int atomicNumberElem // Atomic number of elemental material
138  )
139 {
140  G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
141 
142  G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
143 
144  return ( iter != dedxMapElements.end()) ?
145  (iter->second)->Value( kinEnergyPerNucleon) : 0.0;
146 }
147 
148 // #########################################################################
149 
151  G4double kinEnergyPerNucleon, // Kinetic energy per nucleon
152  G4int atomicNumberIon, // Atomic number of ion
153  const G4String& matIdentifier // Name or chemical formula of material
154  )
155 {
156  G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
157 
158  G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
159 
160  return (iter != dedxMapMaterials.end()) ?
161  (iter->second)->Value(kinEnergyPerNucleon) : 0.0;
162 }
163 
164 // #########################################################################
165 
167  G4PhysicsVector* physicsVector, // Physics vector
168  G4int atomicNumberIon, // Atomic number of ion
169  const G4String& matIdentifier // Name of elemental material
170  )
171 {
172  if(physicsVector == nullptr) {
173  G4Exception ("G4IonStoppingData::AddPhysicsVector() for material",
174  "mat037", FatalException,
175  "Pointer to vector is null-pointer.");
176  return false;
177  }
178 
179  if(matIdentifier.empty()) {
180  G4Exception ("G4IonStoppingData::AddPhysicsVector() for material",
181  "mat038", FatalException, "Invalid name of the material.");
182  return false;
183  }
184 
185  if(atomicNumberIon <= 0) {
186  G4Exception ("G4IonStoppingData::AddPhysicsVector() for material",
187  "mat039", FatalException, "Illegal atomic number.");
188  return false;
189  }
190 
191  G4IonDEDXKeyMat mkey = std::make_pair(atomicNumberIon, matIdentifier);
192 
193  if(dedxMapMaterials.count(mkey) == 1) {
195  ed << "Vector with Z1 = " << atomicNumberIon << ", mat = "
196  << matIdentifier
197  << "already exists. Remove first before replacing.";
198  G4Exception ("G4IonStoppingData::AddPhysicsVector() for material",
199  "mat040", FatalException, ed);
200  return false;
201  }
202 
203  dedxMapMaterials[mkey] = physicsVector;
204 
205  return true;
206 }
207 
208 // #########################################################################
209 
211  G4PhysicsVector* physicsVector, // Physics vector
212  G4int atomicNumberIon, // Atomic number of ion
213  G4int atomicNumberElem // Atomic number of elemental material
214  )
215 {
216  if(physicsVector == nullptr) {
217  G4Exception ("G4IonStoppingData::AddPhysicsVector() for element", "mat037",
218  FatalException, "Pointer to vector is null-pointer.");
219  return false;
220  }
221 
222  if(atomicNumberIon <= 0) {
223  G4Exception ("G4IonStoppingData::AddPhysicsVector() for element", "mat038",
224  FatalException, "Invalid ion number.");
225  return false;
226  }
227 
228  if(atomicNumberElem <= 0) {
229  G4Exception ("G4IonStoppingData::AddPhysicsVector() for element", "mat039",
230  FatalException, "Illegal atomic number.");
231  return false;
232  }
233 
234  G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
235 
236  if(dedxMapElements.count(key) == 1) {
238  ed << "Vector with Z1 = " << atomicNumberIon << ", Z= "
239  << atomicNumberElem
240  << "already exists. Remove first before replacing.";
241  G4Exception ("G4IonStoppingData::AddPhysicsVector() for element", "mat040",
242  FatalException, ed);
243  return false;
244  }
245 
246  dedxMapElements[key] = physicsVector;
247 
248  return true;
249 }
250 
251 // #########################################################################
252 
254  G4int atomicNumberIon, // Atomic number of ion
255  const G4String& matIdentifier // Name or chemical formula of material
256  ) {
257 
258  G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
259 
260  G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
261 
262  if(iter == dedxMapMaterials.end()) {
263  G4Exception ("G4IonStoppingData::RemovePhysicsVector() for material",
264  "mat038", FatalException, "Invalid name of the material.");
265  return false;
266  }
267 
268  G4PhysicsVector* physicsVector = (*iter).second;
269 
270  // Deleting key of physics vector from material map
271  dedxMapMaterials.erase(key);
272 
273  // Deleting physics vector
274  delete physicsVector;
275 
276  return true;
277 }
278 
279 // #########################################################################
280 
282  G4int atomicNumberIon, // Atomic number of ion
283  G4int atomicNumberElem // Atomic number of elemental material
284  ) {
285  G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
286 
287  G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
288 
289  if(iter == dedxMapElements.end()) {
290  G4Exception ("G4IonStoppingData::RemovePhysicsVector() for element",
291  "mat038", FatalException, "Invalid element.");
292  return false;
293  }
294 
295  G4PhysicsVector* physicsVector = (*iter).second;
296 
297  // Deleting key of physics vector from material map
298  dedxMapElements.erase(key);
299 
300  // Deleting physics vector
301  delete physicsVector;
302 
303  return true;
304 }
305 
306 // #########################################################################
307 
309  G4int atomicNumberIon, // Atomic number of ion
310  const G4String& matIdentifier // Name of material
311  ) {
312 
313  if( IsApplicable(atomicNumberIon, matIdentifier) ) return true;
314 
315  char* path = getenv("G4LEDATA");
316  if ( !path ) {
317  G4Exception("G4IonStoppingData::BuildPhysicsVector()", "mat521",
318  FatalException, "G4LEDATA environment variable not set");
319  return false;
320  }
321 
322  std::ostringstream file;
323 
324  file << path << "/" << subDir << "/z"
325  << atomicNumberIon << "_" << matIdentifier
326  << ".dat";
327 
328  G4String fileName = G4String( file.str().c_str() );
329 
330  std::ifstream ifilestream( fileName );
331 
332  if ( !ifilestream.is_open() ) return false;
333 
334  G4LPhysicsFreeVector* physicsVector = new G4LPhysicsFreeVector();
335 
336  if( !physicsVector -> Retrieve(ifilestream, true) ) {
337 
338  ifilestream.close();
339  return false;
340  }
341 
342  physicsVector -> ScaleVector( MeV, MeV * cm2 /( 0.001 * g) );
343  physicsVector -> SetSpline( true );
344  physicsVector -> FillSecondDerivatives();
345 
346  // Retrieved vector is added to material store
347  if( !AddPhysicsVector(physicsVector, atomicNumberIon, matIdentifier) ) {
348  delete physicsVector;
349  ifilestream.close();
350  return false;
351  }
352 
353  ifilestream.close();
354  return true;
355 }
356 
357 // #########################################################################
358 
360  G4int atomicNumberIon, // Atomic number of ion
361  G4int atomicNumberElem // Atomic number of elemental material
362  )
363 {
364  if( IsApplicable(atomicNumberIon, atomicNumberElem) ) return true;
365 
366  char* path = getenv("G4LEDATA");
367  if ( !path ) {
368  G4Exception("G4IonStoppingData::BuildPhysicsVector()", "mat522",
369  FatalException, "G4LEDATA environment variable not set");
370  return false;
371  }
372  std::ostringstream file;
373 
374  file << path << "/" << subDir << "/z"
375  << atomicNumberIon << "_" << atomicNumberElem
376  << ".dat";
377 
378  G4String fileName = G4String( file.str().c_str() );
379 
380  std::ifstream ifilestream( fileName );
381 
382  if ( !ifilestream.is_open() ) return false;
383 
384  G4LPhysicsFreeVector* physicsVector = new G4LPhysicsFreeVector();
385 
386  if( !physicsVector -> Retrieve(ifilestream, true) ) {
387 
388  ifilestream.close();
389  return false;
390  }
391 
392  physicsVector -> ScaleVector( MeV, MeV * cm2 /( 0.001 * g) );
393  physicsVector -> SetSpline( true );
394  physicsVector -> FillSecondDerivatives();
395 
396  // Retrieved vector is added to material store
397  if( !AddPhysicsVector(physicsVector, atomicNumberIon, atomicNumberElem) ) {
398  delete physicsVector;
399  ifilestream.close();
400  return false;
401  }
402 
403  ifilestream.close();
404  return true;
405 }
406 
407 // #########################################################################
408 
410 
411  G4IonDEDXMapMat::iterator iterMat = dedxMapMaterials.begin();
412  G4IonDEDXMapMat::iterator iterMat_end = dedxMapMaterials.end();
413 
414  for(;iterMat != iterMat_end; iterMat++) {
415 
416  G4PhysicsVector* vec = iterMat -> second;
417 
418  if(vec != 0) delete vec;
419  }
420 
421  dedxMapMaterials.clear();
422 
423  G4IonDEDXMapElem::iterator iterElem = dedxMapElements.begin();
424  G4IonDEDXMapElem::iterator iterElem_end = dedxMapElements.end();
425 
426  for(;iterElem != iterElem_end; iterElem++) {
427 
428  G4PhysicsVector* vec = iterElem -> second;
429 
430  if(vec != 0) delete vec;
431  }
432 
433  dedxMapElements.clear();
434 }
435 
436 // #########################################################################
437 
439 
440  G4IonDEDXMapMat::iterator iterMat = dedxMapMaterials.begin();
441  G4IonDEDXMapMat::iterator iterMat_end = dedxMapMaterials.end();
442 
443  G4cout << std::setw(15) << std::right
444  << "Atomic nmb ion"
445  << std::setw(25) << std::right
446  << "Material name"
447  << G4endl;
448 
449  for(;iterMat != iterMat_end; iterMat++) {
450  G4IonDEDXKeyMat key = iterMat -> first;
451  G4PhysicsVector* physicsVector = iterMat -> second;
452 
453  G4int atomicNumberIon = key.first;
454  G4String matIdentifier = key.second;
455 
456  if(physicsVector != 0) {
457  G4cout << std::setw(15) << std::right
458  << atomicNumberIon
459  << std::setw(25) << std::right
460  << matIdentifier
461  << G4endl;
462  }
463  }
464 
465  G4IonDEDXMapElem::iterator iterElem = dedxMapElements.begin();
466  G4IonDEDXMapElem::iterator iterElem_end = dedxMapElements.end();
467 
468  G4cout << std::setw(15) << std::right
469  << "Atomic nmb ion"
470  << std::setw(25) << std::right
471  << "Atomic nmb material"
472  << G4endl;
473 
474  for(;iterElem != iterElem_end; iterElem++) {
475  G4IonDEDXKeyElem key = iterElem -> first;
476  G4PhysicsVector* physicsVector = iterElem -> second;
477 
478  G4int atomicNumberIon = key.first;
479  G4int atomicNumberElem = key.second;
480 
481  if(physicsVector != 0) {
482  G4cout << std::setw(15) << std::right
483  << atomicNumberIon
484  << std::setw(25) << std::right
485  << atomicNumberElem
486  << G4endl;
487  }
488  }
489 
490 }
491 
492 // #########################################################################
493 
G4bool BuildPhysicsVector(G4int ionZ, const G4String &matName)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
static constexpr double cm2
Definition: G4SIunits.hh:120
G4IonDEDXMapMat dedxMapMaterials
static constexpr double second
Definition: G4SIunits.hh:157
static constexpr double g
Definition: G4SIunits.hh:183
int G4int
Definition: G4Types.hh:78
G4bool AddPhysicsVector(G4PhysicsVector *physicsVector, G4int atomicNumberIon, const G4String &matIdentifier)
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
G4IonStoppingData(const G4String &leDirectory)
G4IonDEDXMapElem dedxMapElements
G4bool IsApplicable(G4int atomicNumberIon, G4int atomicNumberElem)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4bool RemovePhysicsVector(G4int atomicNumberIon, const G4String &matIdentifier)
virtual ~G4IonStoppingData()
G4double GetDEDX(G4double kinEnergyPerNucleon, G4int atomicNumberIon, G4int atomicNumberElem)
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
std::pair< G4int, G4int > G4IonDEDXKeyElem
G4PhysicsVector * GetPhysicsVector(G4int atomicNumberIon, G4int atomicNumberElem)
std::pair< G4int, G4String > G4IonDEDXKeyMat