source: trunk/source/materials/src/G4IonStoppingData.cc@ 1325

Last change on this file since 1325 was 1197, checked in by garnier, 16 years ago

nvx fichiers dans CVS

File size: 15.6 KB
Line 
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//
27// ===========================================================================
28// GEANT4 class source file
29//
30// Class: G4IonStoppingData
31//
32// Base class: G4VIonDEDXTable
33//
34// Author: Anton Lechner (Anton.Lechner@cern.ch)
35//
36// First implementation: 03. 11. 2009
37//
38// Modifications:
39//
40//
41// Class description: Class which can read ion stopping power data from
42// $G4LEDATA/ion_stopping_data
43//
44// Comments:
45//
46// ===========================================================================
47//
48
49#include "G4IonStoppingData.hh"
50#include "G4PhysicsVector.hh"
51#include "G4LPhysicsFreeVector.hh"
52#include <fstream>
53#include <sstream>
54#include <iomanip>
55
56
57// #########################################################################
58
59G4IonStoppingData::G4IonStoppingData(const G4String& leDirectory) :
60 subDir( leDirectory ) {
61
62}
63
64// #########################################################################
65
66G4IonStoppingData::~G4IonStoppingData() {
67
68 ClearTable();
69}
70
71// #########################################################################
72
73G4bool G4IonStoppingData::IsApplicable(
74 G4int atomicNumberIon, // Atomic number of ion
75 G4int atomicNumberElem // Atomic number of elemental material
76 ) {
77 G4bool isApplicable = true;
78 G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
79
80 G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
81
82 if(iter == dedxMapElements.end()) isApplicable = false;
83
84 return isApplicable;
85}
86
87// #########################################################################
88
89G4bool G4IonStoppingData::IsApplicable(
90 G4int atomicNumberIon, // Atomic number of ion
91 const G4String& matIdentifier // Name or chemical formula of material
92 ) {
93 G4bool isApplicable = true;
94 G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
95
96 G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
97
98 if(iter == dedxMapMaterials.end()) isApplicable = false;
99
100 return isApplicable;
101}
102
103// #########################################################################
104
105G4PhysicsVector* G4IonStoppingData::GetPhysicsVector(
106 G4int atomicNumberIon, // Atomic number of ion
107 G4int atomicNumberElem // Atomic number of elemental material
108 ) {
109
110 G4PhysicsVector* physVector = 0;
111
112 G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
113
114 G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
115
116 if(iter != dedxMapElements.end()) physVector = iter -> second;
117
118 return physVector;
119}
120
121// #########################################################################
122
123G4PhysicsVector* G4IonStoppingData::GetPhysicsVector(
124 G4int atomicNumberIon, // Atomic number of ion
125 const G4String& matIdentifier // Name or chemical formula of material
126 ) {
127
128 G4PhysicsVector* physVector = 0;
129
130 G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
131
132 G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
133
134 if(iter != dedxMapMaterials.end()) physVector = iter -> second;
135
136 return physVector;
137}
138
139// #########################################################################
140
141G4double G4IonStoppingData::GetDEDX(
142 G4double kinEnergyPerNucleon, // Kinetic energy per nucleon
143 G4int atomicNumberIon, // Atomic number of ion
144 G4int atomicNumberElem // Atomic number of elemental material
145 ) {
146 G4double dedx = 0;
147
148 G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
149
150 G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
151
152 if( iter != dedxMapElements.end() ) {
153 G4PhysicsVector* physicsVector = iter -> second;
154
155 G4bool b;
156 dedx = physicsVector -> GetValue( kinEnergyPerNucleon, b );
157 }
158
159 return dedx;
160}
161
162// #########################################################################
163
164G4double G4IonStoppingData::GetDEDX(
165 G4double kinEnergyPerNucleon, // Kinetic energy per nucleon
166 G4int atomicNumberIon, // Atomic number of ion
167 const G4String& matIdentifier // Name or chemical formula of material
168 ) {
169 G4double dedx = 0;
170
171 G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
172
173 G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
174
175 if(iter != dedxMapMaterials.end()) {
176 G4PhysicsVector* physicsVector = iter -> second;
177
178 G4bool b;
179 dedx = physicsVector -> GetValue( kinEnergyPerNucleon, b );
180 }
181
182 return dedx;
183}
184
185// #########################################################################
186
187G4bool G4IonStoppingData::AddPhysicsVector(
188 G4PhysicsVector* physicsVector, // Physics vector
189 G4int atomicNumberIon, // Atomic number of ion
190 const G4String& matIdentifier // Name of elemental material
191 ) {
192
193 if(physicsVector == 0) {
194
195#ifdef G4VERBOSE
196 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: Pointer to vector"
197 << " is null-pointer."
198 << G4endl;
199#endif
200
201 return false;
202 }
203
204 if(matIdentifier.empty()) {
205
206#ifdef G4VERBOSE
207 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
208 << "Cannot add physics vector. Invalid name."
209 << G4endl;
210#endif
211
212 return false;
213 }
214
215 if(atomicNumberIon <= 0) {
216
217#ifdef G4VERBOSE
218 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
219 << "Cannot add physics vector. Illegal atomic number."
220 << G4endl;
221#endif
222
223 return false;
224 }
225
226 G4IonDEDXKeyMat mkey = std::make_pair(atomicNumberIon, matIdentifier);
227
228 if(dedxMapMaterials.count(mkey) == 1) {
229
230#ifdef G4VERBOSE
231 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
232 << "Vector with Z1 = " << atomicNumberIon << ", mat = "
233 << matIdentifier
234 << "already exists. Remove first before replacing."
235 << G4endl;
236#endif
237
238 return false;
239 }
240
241 dedxMapMaterials[mkey] = physicsVector;
242
243 return true;
244}
245
246// #########################################################################
247
248G4bool G4IonStoppingData::AddPhysicsVector(
249 G4PhysicsVector* physicsVector, // Physics vector
250 G4int atomicNumberIon, // Atomic number of ion
251 G4int atomicNumberElem // Atomic number of elemental material
252 ) {
253
254 if(physicsVector == 0) {
255
256#ifdef G4VERBOSE
257 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
258 << "Pointer to vector is null-pointer."
259 << G4endl;
260#endif
261
262 return false;
263 }
264
265 if(atomicNumberIon <= 0) {
266
267#ifdef G4VERBOSE
268 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
269 << "Cannot add physics vector. Illegal atomic number."
270 << G4endl;
271#endif
272
273 return false;
274 }
275
276 if(atomicNumberElem <= 0) {
277
278#ifdef G4VERBOSE
279 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
280 << "Atomic number of element < 0."
281 << G4endl;
282#endif
283 return false;
284 }
285
286 G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
287
288 if(dedxMapElements.count(key) == 1) {
289
290#ifdef G4VERBOSE
291 G4cerr << "G4IonStoppingData::AddPhysicsVector() Error: "
292 << "Vector with Z1 = " << atomicNumberIon << ", Z2 = "
293 << atomicNumberElem
294 << " already exists. Remove first before replacing."
295 << G4endl;
296#endif
297 return false;
298 }
299
300 dedxMapElements[key] = physicsVector;
301
302 return true;
303}
304
305// #########################################################################
306
307G4bool G4IonStoppingData::RemovePhysicsVector(
308 G4int atomicNumberIon, // Atomic number of ion
309 const G4String& matIdentifier // Name or chemical formula of material
310 ) {
311
312 G4IonDEDXKeyMat key = std::make_pair(atomicNumberIon, matIdentifier);
313
314 G4IonDEDXMapMat::iterator iter = dedxMapMaterials.find(key);
315
316 if(iter == dedxMapMaterials.end()) {
317
318#ifdef G4VERBOSE
319 G4cerr << "G4IonStoppingData::RemovePhysicsVector() Warning: "
320 << "Cannot remove physics vector. Vector not found."
321 << G4endl;
322#endif
323
324 return false;
325 }
326
327 G4PhysicsVector* physicsVector = (*iter).second;
328
329 // Deleting key of physics vector from material map
330 dedxMapMaterials.erase(key);
331
332 // Deleting physics vector
333 delete physicsVector;
334
335 return true;
336}
337
338// #########################################################################
339
340G4bool G4IonStoppingData::RemovePhysicsVector(
341 G4int atomicNumberIon, // Atomic number of ion
342 G4int atomicNumberElem // Atomic number of elemental material
343 ) {
344 G4IonDEDXKeyElem key = std::make_pair(atomicNumberIon, atomicNumberElem);
345
346 G4IonDEDXMapElem::iterator iter = dedxMapElements.find(key);
347
348 if(iter == dedxMapElements.end()) {
349
350#ifdef G4VERBOSE
351 G4cerr << "G4IonStoppingData::RemovePhysicsVector() Warning: "
352 << "Cannot remove physics vector. Vector not found."
353 << G4endl;
354#endif
355
356 return false;
357 }
358
359 G4PhysicsVector* physicsVector = (*iter).second;
360
361 // Deleting key of physics vector from material map
362 dedxMapElements.erase(key);
363
364 // Deleting physics vector
365 delete physicsVector;
366
367 return true;
368}
369
370// #########################################################################
371
372G4bool G4IonStoppingData::BuildPhysicsVector(
373 G4int atomicNumberIon, // Atomic number of ion
374 const G4String& matIdentifier // Name of material
375 ) {
376
377 if( IsApplicable(atomicNumberIon, matIdentifier) ) return true;
378
379 char* path = getenv("G4LEDATA");
380 if ( !path ) G4Exception("G4LEDATA environment variable not set");
381
382 std::ostringstream file;
383
384 file << path << "/" << subDir << "/z"
385 << atomicNumberIon << "_" << matIdentifier
386 << ".dat";
387
388 G4String fileName = G4String( file.str().c_str() );
389
390 std::ifstream ifilestream( fileName );
391
392 if ( !ifilestream.is_open() ) return false;
393
394 G4LPhysicsFreeVector* physicsVector = new G4LPhysicsFreeVector();
395
396 if( !physicsVector -> Retrieve(ifilestream, true) ) {
397
398 ifilestream.close();
399 return false;
400 }
401
402 physicsVector -> ScaleVector( MeV, MeV * cm2 /( 0.001 * g) );
403 physicsVector -> SetSpline( true );
404 physicsVector -> FillSecondDerivatives();
405
406 // Retrieved vector is added to material store
407 if( !AddPhysicsVector(physicsVector, atomicNumberIon, matIdentifier) ) {
408 delete physicsVector;
409 ifilestream.close();
410 return false;
411 }
412
413 ifilestream.close();
414 return true;
415}
416
417// #########################################################################
418
419G4bool G4IonStoppingData::BuildPhysicsVector(
420 G4int atomicNumberIon, // Atomic number of ion
421 G4int atomicNumberElem // Atomic number of elemental material
422 ) {
423
424 if( IsApplicable(atomicNumberIon, atomicNumberElem) ) return true;
425
426 char* path = getenv("G4LEDATA");
427 if ( !path ) G4Exception("G4LEDATA environment variable not set");
428
429 std::ostringstream file;
430
431 file << path << "/" << subDir << "/z"
432 << atomicNumberIon << "_" << atomicNumberElem
433 << ".dat";
434
435 G4String fileName = G4String( file.str().c_str() );
436
437 std::ifstream ifilestream( fileName );
438
439 if ( !ifilestream.is_open() ) return false;
440
441 G4LPhysicsFreeVector* physicsVector = new G4LPhysicsFreeVector();
442
443 if( !physicsVector -> Retrieve(ifilestream, true) ) {
444
445 ifilestream.close();
446 return false;
447 }
448
449 physicsVector -> ScaleVector( MeV, MeV * cm2 /( 0.001 * g) );
450 physicsVector -> SetSpline( true );
451 physicsVector -> FillSecondDerivatives();
452
453 // Retrieved vector is added to material store
454 if( !AddPhysicsVector(physicsVector, atomicNumberIon, atomicNumberElem) ) {
455 delete physicsVector;
456 ifilestream.close();
457 return false;
458 }
459
460 ifilestream.close();
461 return true;
462}
463
464// #########################################################################
465
466void G4IonStoppingData::ClearTable() {
467
468 G4IonDEDXMapMat::iterator iterMat = dedxMapMaterials.begin();
469 G4IonDEDXMapMat::iterator iterMat_end = dedxMapMaterials.end();
470
471 for(;iterMat != iterMat_end; iterMat++) {
472
473 G4PhysicsVector* vec = iterMat -> second;
474
475 if(vec != 0) delete vec;
476 }
477
478 dedxMapMaterials.clear();
479
480 G4IonDEDXMapElem::iterator iterElem = dedxMapElements.begin();
481 G4IonDEDXMapElem::iterator iterElem_end = dedxMapElements.end();
482
483 for(;iterElem != iterElem_end; iterElem++) {
484
485 G4PhysicsVector* vec = iterElem -> second;
486
487 if(vec != 0) delete vec;
488 }
489
490 dedxMapElements.clear();
491}
492
493// #########################################################################
494
495void G4IonStoppingData::DumpMap() {
496
497 G4IonDEDXMapMat::iterator iterMat = dedxMapMaterials.begin();
498 G4IonDEDXMapMat::iterator iterMat_end = dedxMapMaterials.end();
499
500 G4cout << std::setw(15) << std::right
501 << "Atomic nmb ion"
502 << std::setw(25) << std::right
503 << "Material name"
504 << G4endl;
505
506 for(;iterMat != iterMat_end; iterMat++) {
507 G4IonDEDXKeyMat key = iterMat -> first;
508 G4PhysicsVector* physicsVector = iterMat -> second;
509
510 G4int atomicNumberIon = key.first;
511 G4String matIdentifier = key.second;
512
513 if(physicsVector != 0) {
514 G4cout << std::setw(15) << std::right
515 << atomicNumberIon
516 << std::setw(25) << std::right
517 << matIdentifier
518 << G4endl;
519 }
520 }
521
522 G4IonDEDXMapElem::iterator iterElem = dedxMapElements.begin();
523 G4IonDEDXMapElem::iterator iterElem_end = dedxMapElements.end();
524
525 G4cout << std::setw(15) << std::right
526 << "Atomic nmb ion"
527 << std::setw(25) << std::right
528 << "Atomic nmb material"
529 << G4endl;
530
531 for(;iterElem != iterElem_end; iterElem++) {
532 G4IonDEDXKeyElem key = iterElem -> first;
533 G4PhysicsVector* physicsVector = iterElem -> second;
534
535 G4int atomicNumberIon = key.first;
536 G4int atomicNumberElem = key.second;
537
538 if(physicsVector != 0) {
539 G4cout << std::setw(15) << std::right
540 << atomicNumberIon
541 << std::setw(25) << std::right
542 << atomicNumberElem
543 << G4endl;
544 }
545 }
546
547}
548
549// #########################################################################
550
Note: See TracBrowser for help on using the repository browser.