Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DNAMolecularMaterial.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: G4DNAMolecularMaterial.cc 101354 2016-11-15 08:27:51Z gcosmo $
27 //
28 // Author: Mathieu Karamitros
29 //
31 #include "G4Material.hh"
32 #include <utility>
33 #include "G4StateManager.hh"
34 #include "G4Threading.hh"
35 #include "G4AutoLock.hh"
36 #include "G4StateManager.hh"
37 #include "G4MoleculeTable.hh"
38 
39 using namespace std;
40 
43 
44 //------------------------------------------------------------------------------
45 
47  const G4Material* mat2) const
48 {
49  if (mat1 == 0 && mat2 == 0) return false; //(mat1 == mat2)
50  if (mat1 == 0) return true; // mat1 < mat2
51  if (mat2 == 0) return false; //mat2 < mat1
52 
53  const G4Material* baseMat1 = mat1->GetBaseMaterial();
54  const G4Material* baseMat2 = mat2->GetBaseMaterial();
55 
56  if ((baseMat1 || baseMat2) == 0){
57  // None of the materials derives from a base material
58  return mat1 < mat2;
59  }
60  else if (baseMat1 && baseMat2){
61  // Both materials derive from a base material
62  return baseMat1 < baseMat2;
63  }
64 
65  else if (baseMat1 && (baseMat2 == 0)){
66  // Only the material 1 derives from a base material
67  return baseMat1 < mat2;
68  }
69  // only case baseMat1==0 && baseMat2 remains
70  return mat1 < baseMat2;
71 }
72 
73 //------------------------------------------------------------------------------
74 
76 {
77  if (!fInstance) new G4DNAMolecularMaterial();
78  return fInstance;
79 }
80 
81 //------------------------------------------------------------------------------
82 
84 {
85  if (fInstance){
86  delete fInstance;
87  fInstance = 0;
88  }
89 }
90 
91 //------------------------------------------------------------------------------
92 
94 {
95  fpCompFractionTable = 0;
96  fpCompDensityTable = 0;
97  fpCompNumMolPerVolTable = 0;
98  fIsInitialized = false;
99  fNMaterials = 0;
100  fInstance = this;
101 }
102 
103 //------------------------------------------------------------------------------
104 
106 {
107  if (fpCompFractionTable){
108  fpCompFractionTable->clear();
109  delete fpCompFractionTable;
110  fpCompFractionTable = 0;
111  }
112  if (fpCompDensityTable){
113  fpCompDensityTable->clear();
114  delete fpCompDensityTable;
115  fpCompDensityTable = 0;
116  }
117  if (fpCompNumMolPerVolTable){
118  fpCompNumMolPerVolTable->clear();
119  delete fpCompNumMolPerVolTable;
120  fpCompNumMolPerVolTable = 0;
121  }
122 
123  map<const G4Material*, std::vector<double>*, CompareMaterial>::iterator it;
124 
125  for (it = fAskedDensityTable.begin(); it != fAskedDensityTable.end(); it++){
126  if (it->second){
127  delete it->second;
128  it->second = 0;
129  }
130  }
131 
132  for (it = fAskedNumPerVolTable.begin(); it != fAskedNumPerVolTable.end();
133  it++){
134  if (it->second){
135  delete it->second;
136  it->second = 0;
137  }
138  }
139 }
140 
141 //------------------------------------------------------------------------------
142 
145 {
146  Create();
147  fInstance = this;
148 }
149 
150 //------------------------------------------------------------------------------
151 
153 {
154  if (requestedState == G4State_Idle && G4StateManager::GetStateManager()
155  ->GetPreviousState() == G4State_PreInit){
156  Initialize();
157  }
158  else if (requestedState == G4State_Quit){
159 // G4cout << "G4DNAMolecularMaterial::Notify ---> received G4State_Quit"
160 // << G4endl;
161  Clear();
162  //DeleteInstance();
163  }
164  return true;
165 }
166 
167 //------------------------------------------------------------------------------
168 
170  const G4DNAMolecularMaterial& /*rhs*/) :
172 {
173  Create();
174 }
175 
176 //------------------------------------------------------------------------------
177 
180 {
181  if (this == &rhs) return *this;
182  Create();
183  return *this;
184 }
185 
186 //------------------------------------------------------------------------------
187 
189 {
190 // G4cout << "Deleting G4DNAMolecularMaterial" << G4endl;
191  Clear();
192  fInstance = 0;
193  //assert(G4StateManager::GetStateManager()->DeregisterDependent(this) == true);
194 }
195 
196 //------------------------------------------------------------------------------
197 
199 {
200  G4AutoLock l(&aMutex);
201  if (fIsInitialized){
202  return;
203  }
204 
205  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
206 
207  fNMaterials = materialTable->size();
208  // This is to prevent segment fault if materials are created later on
209  // Actually this creation should not be done
210 
211  if (fpCompFractionTable == 0){
212  fpCompFractionTable = new vector<ComponentMap>(materialTable->size());
213  }
214 
215  G4Material* mat(0);
216 
217  for (size_t i = 0; i < fNMaterials; i++){
218  mat = materialTable->at(i);
219  SearchMolecularMaterial(mat, mat, 1);
220 
221  mat = 0;
222  }
223 
226 
227  fIsInitialized = true;
228 }
229 
230 //------------------------------------------------------------------------------
231 
233 {
234  if (fpCompFractionTable){
235  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
236  fpCompDensityTable = new vector<ComponentMap>(
237  G4Material::GetMaterialTable()->size());
238 
239  G4Material* parentMat;
240  const G4Material* compMat(0);
241  double massFraction = -1;
242  double parentDensity = -1;
243 
244  for (size_t i = 0; i < fNMaterials; i++){
245  parentMat = materialTable->at(i);
246  ComponentMap& massFractionComp = (*fpCompFractionTable)[i];
247  ComponentMap& densityComp = (*fpCompDensityTable)[i];
248 
249  parentDensity = parentMat->GetDensity();
250 
251  for (ComponentMap::iterator it = massFractionComp.begin();
252  it != massFractionComp.end(); it++){
253  compMat = it->first;
254  massFraction = it->second;
255  densityComp[compMat] = massFraction * parentDensity;
256  compMat = 0;
257  massFraction = -1;
258  }
259  }
260  }
261  else{
262  G4ExceptionDescription exceptionDescription;
263  exceptionDescription << "The pointer fpCompFractionTable is not initialized"
264  << G4endl;
265  G4Exception("G4DNAMolecularMaterial::InitializeDensity",
266  "G4DNAMolecularMaterial001", FatalException,
267  exceptionDescription);
268  }
269 }
270 
271 //------------------------------------------------------------------------------
272 
274 {
275  if (fpCompDensityTable){
276  fpCompNumMolPerVolTable = new vector<ComponentMap>(fNMaterials);
277 
278  const G4Material* compMat(0);
279 
280  for (size_t i = 0; i < fNMaterials; i++){
281  ComponentMap& massFractionComp = (*fpCompFractionTable)[i];
282  ComponentMap& densityComp = (*fpCompDensityTable)[i];
283  ComponentMap& numMolPerVol = (*fpCompNumMolPerVolTable)[i];
284 
285  for (ComponentMap::iterator it = massFractionComp.begin();
286  it != massFractionComp.end(); it++){
287  compMat = it->first;
288  numMolPerVol[compMat] = densityComp[compMat]
289  / compMat->GetMassOfMolecule();
290  compMat = 0;
291  }
292  }
293  }
294  else{
295  G4ExceptionDescription exceptionDescription;
296  exceptionDescription << "The pointer fpCompDensityTable is not initialized"
297  << G4endl;
298  G4Exception("G4DNAMolecularMaterial::InitializeNumMolPerVol",
299  "G4DNAMolecularMaterial002", FatalException,
300  exceptionDescription);
301  }
302 }
303 
304 //------------------------------------------------------------------------------
305 
306 void
308  G4Material* molecularMaterial,
309  G4double fraction)
310 {
311  ComponentMap& matComponent =
312  (*fpCompFractionTable)[parentMaterial->GetIndex()];
313 
314  if (matComponent.empty()){
315  matComponent[molecularMaterial] = fraction;
316  return;
317  }
318 
319  ComponentMap::iterator it = matComponent.find(molecularMaterial);
320 
321  if (it == matComponent.end()){
322  matComponent[molecularMaterial] = fraction;
323  }
324  else{
325  matComponent[molecularMaterial] = it->second + fraction;
326  }
327 }
328 
329 //------------------------------------------------------------------------------
330 
333  double currentFraction)
334 {
335  if (material->GetMassOfMolecule() != 0.0){
336  RecordMolecularMaterial(parentMaterial, material, currentFraction);
337  return;
338  }
339 
340  G4Material* compMat(0);
341  G4double fraction = -1;
342  std::map<G4Material*, G4double> matComponent = material->GetMatComponents();
343  std::map<G4Material*, G4double>::iterator it = matComponent.begin();
344 
345  for (; it != matComponent.end(); it++){
346  compMat = it->first;
347  fraction = it->second;
348  if (compMat->GetMassOfMolecule() == 0.0){
349  SearchMolecularMaterial(parentMaterial, compMat,
350  currentFraction * fraction);
351  }
352  else{
353  RecordMolecularMaterial(parentMaterial, compMat,
354  currentFraction * fraction);
355  }
356 
357  //compMat = 0;
358  //fraction = -1;
359  }
360 }
361 
362 //------------------------------------------------------------------------------
363 
364 const std::vector<double>*
366 GetDensityTableFor(const G4Material* lookForMaterial) const
367 {
368  if (!fpCompDensityTable){
369  if (fIsInitialized){
370  G4ExceptionDescription exceptionDescription;
371  exceptionDescription
372  << "The pointer fpCompDensityTable is not initialized will the "
373  "singleton of G4DNAMolecularMaterial "
374  << "has already been initialized." << G4endl;
375  G4Exception("G4DNAMolecularMaterial::GetDensityTableFor",
376  "G4DNAMolecularMaterial003", FatalException,
377  exceptionDescription);
378  }
379 
380  if (G4StateManager::GetStateManager()->GetCurrentState() == G4State_Idle){
381  const_cast<G4DNAMolecularMaterial*>(this)->Initialize();
382  }
383  else{
384  G4ExceptionDescription exceptionDescription;
385  exceptionDescription
386  << "The geant4 application is at the wrong state. State must be: "
387  "G4State_Idle."
388  << G4endl;
389  G4Exception("G4DNAMolecularMaterial::GetDensityTableFor",
390  "G4DNAMolecularMaterial_WRONG_STATE_APPLICATION",
391  FatalException, exceptionDescription);
392  }
393  }
394 
395  std::map<const G4Material*, std::vector<double>*, CompareMaterial>::
396  const_iterator it_askedDensityTable =
397  fAskedDensityTable.find(lookForMaterial);
398 
399  if (it_askedDensityTable != fAskedDensityTable.end()){
400  return it_askedDensityTable->second;
401  }
402 
403  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
404 
405  std::vector<double>* output = new std::vector<double>(materialTable->size());
406 
407  ComponentMap::const_iterator it;
408 
409  G4bool materialWasNotFound = true;
410 
411  for (size_t i = 0; i < fNMaterials; i++){
412  ComponentMap& densityTable = (*fpCompDensityTable)[i];
413 
414  it = densityTable.find(lookForMaterial);
415 
416  if (it == densityTable.end()){
417  (*output)[i] = 0.0;
418  }
419  else{
420  materialWasNotFound = false;
421  (*output)[i] = it->second;
422  }
423  }
424 
425  if (materialWasNotFound){
426  PrintNotAMolecularMaterial("G4DNAMolecularMaterial::GetDensityTableFor",
427  lookForMaterial);
428  }
429 
430  fAskedDensityTable.insert(make_pair(lookForMaterial, output));
431 
432  return output;
433 }
434 
435 //------------------------------------------------------------------------------
436 
438  const G4Material* lookForMaterial) const
439 {
440  if(lookForMaterial==0) return nullptr;
441 
443  if (fIsInitialized){
444  G4ExceptionDescription exceptionDescription;
445  exceptionDescription
446  << "The pointer fpCompNumMolPerVolTable is not initialized whereas "
447  "the singleton of G4DNAMolecularMaterial "
448  << "has already been initialized." << G4endl;
449  G4Exception("G4DNAMolecularMaterial::GetNumMolPerVolTableFor",
450  "G4DNAMolecularMaterial005", FatalException,
451  exceptionDescription);
452  }
453 
454  if (G4StateManager::GetStateManager()->GetCurrentState() == G4State_Idle){
455  const_cast<G4DNAMolecularMaterial*>(this)->Initialize();
456  }
457  else{
458  G4ExceptionDescription exceptionDescription;
459  exceptionDescription
460  << "The geant4 application is at the wrong state. State must be : "
461  "G4State_Idle."
462  << G4endl;
463  G4Exception("G4DNAMolecularMaterial::GetNumMolPerVolTableFor",
464  "G4DNAMolecularMaterial_WRONG_STATE_APPLICATION",
465  FatalException, exceptionDescription);
466  }
467  }
468 
469  std::map<const G4Material*, std::vector<double>*, CompareMaterial>::
470  const_iterator it_askedNumMolPerVolTable =
471  fAskedNumPerVolTable.find(lookForMaterial);
472  if (it_askedNumMolPerVolTable != fAskedNumPerVolTable.end()){
473  return it_askedNumMolPerVolTable->second;
474  }
475 
476  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
477 
478  std::vector<double>* output = new std::vector<double>(materialTable->size());
479 
480  ComponentMap::const_iterator it;
481 
482  G4bool materialWasNotFound = true;
483 
484  for (size_t i = 0; i < fNMaterials; i++){
485  ComponentMap& densityTable = (*fpCompNumMolPerVolTable)[i];
486 
487  it = densityTable.find(lookForMaterial);
488 
489  if (it == densityTable.end()){
490  (*output)[i] = 0.0;
491  }
492  else{
493  materialWasNotFound = false;
494  (*output)[i] = it->second;
495  }
496  }
497 
498  if (materialWasNotFound){
500  "G4DNAMolecularMaterial::GetNumMolPerVolTableFor", lookForMaterial);
501  }
502 
503  fAskedNumPerVolTable.insert(make_pair(lookForMaterial, output));
504 
505  return output;
506 }
507 
508 //------------------------------------------------------------------------------
509 
511 PrintNotAMolecularMaterial(const char* methodName,
512  const G4Material* lookForMaterial) const
513 {
514  std::map<const G4Material*, bool, CompareMaterial>::iterator it =
515  fWarningPrinted.find(lookForMaterial);
516 
517  if (it == fWarningPrinted.end()){
518  G4ExceptionDescription exceptionDescription;
519  exceptionDescription << "The material " << lookForMaterial->GetName()
520  << " is not defined as a molecular material."
521  << G4endl
522  << "Meaning: The elements should be added to the "
523  "material using atom count rather than mass fraction "
524  "(cf. G4Material)"
525  << G4endl
526  << "If you want to use DNA processes on liquid water, you should better use "
527  "the NistManager to create the water material."
528  << G4endl
529  << "Since this message is displayed, it means that the DNA models will not "
530  "be called."
531  << "Please note that this message will only appear once even if you are "
532  "using other methods of G4DNAMolecularMaterial."
533  << G4endl;
534 
535  G4Exception(methodName, "MATERIAL_NOT_DEFINE_USING_ATOM_COUNT", JustWarning,
536  exceptionDescription);
537  fWarningPrinted[lookForMaterial] = true;
538  }
539 }
540 
541 //------------------------------------------------------------------------------
542 
546 {
547  int material_id = material->GetIndex();
548  auto it = fMaterialToMolecularConf.find(material_id);
549  if(it == fMaterialToMolecularConf.end()) return 0;
550  return it->second;
551 }
552 
553 //------------------------------------------------------------------------------
554 
555 void
558  G4MolecularConfiguration* molConf)
559 {
560  assert(material != 0);
561  int material_id = material->GetIndex();
562  fMaterialToMolecularConf[material_id] = molConf;
563 }
564 
565 //------------------------------------------------------------------------------
566 
567 void
569  const G4String& molUserID)
570 {
571  assert(material != 0);
572  int material_id = material->GetIndex();
573  fMaterialToMolecularConf[material_id] =
574  G4MoleculeTable::Instance()->GetConfiguration(molUserID, true);
575 }
576 
577 //------------------------------------------------------------------------------
578 
579 void
581  const G4String& molUserID)
582 {
584 
585  if(material == 0){
586  G4cout<< "Material " << materialName
587  << " was not found and therefore won't be linked to "
588  << molUserID << G4endl;
589  return;
590  }
591  SetMolecularConfiguration(material, molUserID);
592 }
const std::vector< double > * GetDensityTableFor(const G4Material *) const
const std::map< G4Material *, G4double > & GetMatComponents() const
Definition: G4Material.hh:237
std::map< const G4Material *, double, CompareMaterial > ComponentMap
void SearchMolecularMaterial(G4Material *parentMaterial, G4Material *material, double currentFraction)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void SetMolecularConfiguration(const G4Material *, G4MolecularConfiguration *)
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:602
virtual G4bool Notify(G4ApplicationState requestedState)
size_t GetIndex() const
Definition: G4Material.hh:262
std::map< const G4Material *, std::vector< double > *, CompareMaterial > fAskedNumPerVolTable
const G4String & GetName() const
Definition: G4Material.hh:178
void RecordMolecularMaterial(G4Material *parentMaterial, G4Material *molecularMaterial, G4double fraction)
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
std::vector< G4Material * > G4MaterialTable
bool operator()(const G4Material *mat1, const G4Material *mat2) const
void PrintNotAMolecularMaterial(const char *methodName, const G4Material *lookForMaterial) const
std::map< int, G4MolecularConfiguration * > fMaterialToMolecularConf
#define G4MUTEX_INITIALIZER
Definition: G4Threading.hh:175
std::vector< ComponentMap > * fpCompNumMolPerVolTable
std::map< const G4Material *, std::vector< double > *, CompareMaterial > fAskedDensityTable
string material
Definition: eplot.py:19
static G4DNAMolecularMaterial * fInstance
static G4StateManager * GetStateManager()
G4GLOB_DLL std::ostream G4cout
const std::vector< double > * GetNumMolPerVolTableFor(const G4Material *) const
bool G4bool
Definition: G4Types.hh:79
G4DNAMolecularMaterial & operator=(const G4DNAMolecularMaterial &)
G4Mutex aMutex
static G4MoleculeTable * Instance()
std::map< const G4Material *, bool, CompareMaterial > fWarningPrinted
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4int G4Mutex
Definition: G4Threading.hh:173
static G4DNAMolecularMaterial * Instance()
std::vector< ComponentMap > * fpCompFractionTable
G4double GetMassOfMolecule() const
Definition: G4Material.hh:242
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:233
std::vector< ComponentMap > * fpCompDensityTable
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4ApplicationState
G4MolecularConfiguration * GetMolecularConfiguration(const G4Material *) const
G4MolecularConfiguration * GetConfiguration(const G4String &, bool mustExist=true)