ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4GDMLWriteMaterials.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4GDMLWriteMaterials.cc
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 // class G4GDMLWriteMaterials Implementation
29 //
30 // Original author: Zoltan Torzsok, November 2007
31 //
32 // --------------------------------------------------------------------
33 
34 #include <sstream>
35 
36 #include "G4GDMLWriteMaterials.hh"
37 
38 #include "G4PhysicalConstants.hh"
39 #include "G4SystemOfUnits.hh"
40 #include "G4Element.hh"
41 #include "G4Isotope.hh"
42 #include "G4Material.hh"
43 
45  : G4GDMLWriteDefine(), materialsElement(0)
46 {
47 }
48 
50 {
51 }
52 
54 AtomWrite(xercesc::DOMElement* element,const G4double& a)
55 {
56  xercesc::DOMElement* atomElement = NewElement("atom");
57  atomElement->setAttributeNode(NewAttribute("unit","g/mole"));
58  atomElement->setAttributeNode(NewAttribute("value",a*mole/g));
59  element->appendChild(atomElement);
60 }
61 
63 DWrite(xercesc::DOMElement* element,const G4double& d)
64 {
65  xercesc::DOMElement* DElement = NewElement("D");
66  DElement->setAttributeNode(NewAttribute("unit","g/cm3"));
67  DElement->setAttributeNode(NewAttribute("value",d*cm3/g));
68  element->appendChild(DElement);
69 }
70 
72 PWrite(xercesc::DOMElement* element,const G4double& P)
73 {
74  xercesc::DOMElement* PElement = NewElement("P");
75  PElement->setAttributeNode(NewAttribute("unit","pascal"));
76  PElement->setAttributeNode(NewAttribute("value",P/hep_pascal));
77  element->appendChild(PElement);
78 }
79 
81 TWrite(xercesc::DOMElement* element,const G4double& T)
82 {
83  xercesc::DOMElement* TElement = NewElement("T");
84  TElement->setAttributeNode(NewAttribute("unit","K"));
85  TElement->setAttributeNode(NewAttribute("value",T/kelvin));
86  element->appendChild(TElement);
87 }
88 
90 MEEWrite(xercesc::DOMElement* element,const G4double& MEE)
91 {
92  xercesc::DOMElement* PElement = NewElement("MEE");
93  PElement->setAttributeNode(NewAttribute("unit","eV"));
94  PElement->setAttributeNode(NewAttribute("value",MEE/electronvolt));
95  element->appendChild(PElement);
96 }
97 
99 IsotopeWrite(const G4Isotope* const isotopePtr)
100 {
101  const G4String name = GenerateName(isotopePtr->GetName(),isotopePtr);
102 
103  xercesc::DOMElement* isotopeElement = NewElement("isotope");
104  isotopeElement->setAttributeNode(NewAttribute("name",name));
105  isotopeElement->setAttributeNode(NewAttribute("N",isotopePtr->GetN()));
106  isotopeElement->setAttributeNode(NewAttribute("Z",isotopePtr->GetZ()));
107  materialsElement->appendChild(isotopeElement);
108  AtomWrite(isotopeElement,isotopePtr->GetA());
109 }
110 
111 void G4GDMLWriteMaterials::ElementWrite(const G4Element* const elementPtr)
112 {
113  const G4String name = GenerateName(elementPtr->GetName(),elementPtr);
114 
115  xercesc::DOMElement* elementElement = NewElement("element");
116  elementElement->setAttributeNode(NewAttribute("name",name));
117 
118  const size_t NumberOfIsotopes = elementPtr->GetNumberOfIsotopes();
119 
120  if (NumberOfIsotopes>0)
121  {
122  const G4double* RelativeAbundanceVector =
123  elementPtr->GetRelativeAbundanceVector();
124  for (size_t i=0;i<NumberOfIsotopes;i++)
125  {
126  G4String fractionref = GenerateName(elementPtr->GetIsotope(i)->GetName(),
127  elementPtr->GetIsotope(i));
128  xercesc::DOMElement* fractionElement = NewElement("fraction");
129  fractionElement->setAttributeNode(NewAttribute("n",
130  RelativeAbundanceVector[i]));
131  fractionElement->setAttributeNode(NewAttribute("ref",fractionref));
132  elementElement->appendChild(fractionElement);
133  AddIsotope(elementPtr->GetIsotope(i));
134  }
135  }
136  else
137  {
138  elementElement->setAttributeNode(NewAttribute("Z",elementPtr->GetZ()));
139  AtomWrite(elementElement,elementPtr->GetA());
140  }
141 
142  materialsElement->appendChild(elementElement);
143  // Append the element AFTER all the possible components are appended!
144 }
145 
146 void G4GDMLWriteMaterials::MaterialWrite(const G4Material* const materialPtr)
147 {
148  G4String state_str("undefined");
149  const G4State state = materialPtr->GetState();
150  if (state==kStateSolid) { state_str = "solid"; } else
151  if (state==kStateLiquid) { state_str = "liquid"; } else
152  if (state==kStateGas) { state_str = "gas"; }
153 
154  const G4String name = GenerateName(materialPtr->GetName(), materialPtr);
155 
156  xercesc::DOMElement* materialElement = NewElement("material");
157  materialElement->setAttributeNode(NewAttribute("name",name));
158  materialElement->setAttributeNode(NewAttribute("state",state_str));
159 
160  // Write any property attached to the material...
161  //
162  if (materialPtr->GetMaterialPropertiesTable())
163  {
164  PropertyWrite(materialElement, materialPtr);
165  }
166 
167  if (materialPtr->GetTemperature() != STP_Temperature)
168  { TWrite(materialElement,materialPtr->GetTemperature()); }
169  if (materialPtr->GetPressure() != STP_Pressure)
170  { PWrite(materialElement,materialPtr->GetPressure()); }
171 
172  // Write Ionisation potential (mean excitation energy)
173  MEEWrite(materialElement,materialPtr->GetIonisation()->GetMeanExcitationEnergy());
174 
175  DWrite(materialElement,materialPtr->GetDensity());
176 
177  const size_t NumberOfElements = materialPtr->GetNumberOfElements();
178 
179  if ( (NumberOfElements>1)
180  || ( materialPtr->GetElement(0)
181  && materialPtr->GetElement(0)->GetNumberOfIsotopes()>1 ) )
182  {
183  const G4double* MassFractionVector = materialPtr->GetFractionVector();
184 
185  for (size_t i=0;i<NumberOfElements;i++)
186  {
187  const G4String fractionref =
188  GenerateName(materialPtr->GetElement(i)->GetName(),
189  materialPtr->GetElement(i));
190  xercesc::DOMElement* fractionElement = NewElement("fraction");
191  fractionElement->setAttributeNode(NewAttribute("n",
192  MassFractionVector[i]));
193  fractionElement->setAttributeNode(NewAttribute("ref",fractionref));
194  materialElement->appendChild(fractionElement);
195  AddElement(materialPtr->GetElement(i));
196  }
197  }
198  else
199  {
200  materialElement->setAttributeNode(NewAttribute("Z",materialPtr->GetZ()));
201  AtomWrite(materialElement,materialPtr->GetA());
202  }
203 
204  // Append the material AFTER all the possible components are appended!
205  //
206  materialsElement->appendChild(materialElement);
207 }
208 
210  const G4PhysicsOrderedFreeVector* const pvec)
211 {
212  for (size_t i=0; i<propertyList.size(); i++) // Check if property is
213  { // already in the list!
214  if (propertyList[i] == pvec) { return; }
215  }
216  propertyList.push_back(pvec);
217 
218  const G4String matrixref = GenerateName(key, pvec);
219  xercesc::DOMElement* matrixElement = NewElement("matrix");
220  matrixElement->setAttributeNode(NewAttribute("name", matrixref));
221  matrixElement->setAttributeNode(NewAttribute("coldim", "2"));
222  std::ostringstream pvalues;
223  for (size_t i=0; i<pvec->GetVectorLength(); i++)
224  {
225  if (i!=0) { pvalues << " "; }
226  pvalues << pvec->Energy(i) << " " << (*pvec)[i];
227  }
228  matrixElement->setAttributeNode(NewAttribute("values", pvalues.str()));
229 
230  defineElement->appendChild(matrixElement);
231 }
232 
233 void G4GDMLWriteMaterials::PropertyWrite(xercesc::DOMElement* matElement,
234  const G4Material* const mat)
235 {
236  xercesc::DOMElement* propElement;
238 
239  const std::map< G4int, G4PhysicsOrderedFreeVector*,
240  std::less<G4int> >* pmap = ptable->GetPropertyMap();
241  const std::map< G4int, G4double,
242  std::less<G4int> >* cmap = ptable->GetConstPropertyMap();
243  std::map< G4int, G4PhysicsOrderedFreeVector*,
244  std::less<G4int> >::const_iterator mpos;
245  std::map< G4int, G4double,
246  std::less<G4int> >::const_iterator cpos;
247 
248 
249  for (mpos=pmap->begin(); mpos!=pmap->end(); mpos++)
250  {
251  propElement = NewElement("property");
252  propElement->setAttributeNode(NewAttribute("name",
253  ptable->GetMaterialPropertyNames()[mpos->first]));
254  propElement->setAttributeNode(NewAttribute("ref",
255  GenerateName(ptable->GetMaterialPropertyNames()[mpos->first],
256  mpos->second)));
257  if (mpos->second)
258  {
259  PropertyVectorWrite(ptable->GetMaterialPropertyNames()[mpos->first],
260  mpos->second);
261  matElement->appendChild(propElement);
262  }
263  else
264  {
265  G4String warn_message = "Null pointer for material property -"
266  + ptable->GetMaterialPropertyNames()[mpos->first] + "- of material -"
267  + mat->GetName() + "- !";
268  G4Exception("G4GDMLWriteMaterials::PropertyWrite()", "NullPointer",
269  JustWarning, warn_message);
270  continue;
271  }
272  }
273 
274  for (cpos=cmap->begin(); cpos!=cmap->end(); cpos++)
275  {
276  propElement = NewElement("property");
277  propElement->setAttributeNode(NewAttribute("name",
278  ptable->GetMaterialConstPropertyNames()[cpos->first]));
279  propElement->setAttributeNode(NewAttribute("ref",
280  ptable->GetMaterialConstPropertyNames()[cpos->first]));
281  xercesc::DOMElement* constElement = NewElement("constant");
282  constElement->setAttributeNode(NewAttribute("name",
283  ptable->GetMaterialConstPropertyNames()[cpos->first]));
284  constElement->setAttributeNode(NewAttribute("value", cpos->second));
285  defineElement->appendChild(constElement);
286  matElement->appendChild(propElement);
287  }
288 }
289 
290 void G4GDMLWriteMaterials::MaterialsWrite(xercesc::DOMElement* element)
291 {
292 #ifdef G4VERBOSE
293  G4cout << "G4GDML: Writing materials..." << G4endl;
294 #endif
295  materialsElement = NewElement("materials");
296  element->appendChild(materialsElement);
297 
298  isotopeList.clear();
299  elementList.clear();
300  materialList.clear();
301  propertyList.clear();
302 }
303 
304 void G4GDMLWriteMaterials::AddIsotope(const G4Isotope* const isotopePtr)
305 {
306  for (size_t i=0; i<isotopeList.size(); i++) // Check if isotope is
307  { // already in the list!
308  if (isotopeList[i] == isotopePtr) { return; }
309  }
310  isotopeList.push_back(isotopePtr);
311  IsotopeWrite(isotopePtr);
312 }
313 
314 void G4GDMLWriteMaterials::AddElement(const G4Element* const elementPtr)
315 {
316  for (size_t i=0;i<elementList.size();i++) // Check if element is
317  { // already in the list!
318  if (elementList[i] == elementPtr) { return; }
319  }
320  elementList.push_back(elementPtr);
321  ElementWrite(elementPtr);
322 }
323 
324 void G4GDMLWriteMaterials::AddMaterial(const G4Material* const materialPtr)
325 {
326  for (size_t i=0;i<materialList.size();i++) // Check if material is
327  { // already in the list!
328  if (materialList[i] == materialPtr) { return; }
329  }
330  materialList.push_back(materialPtr);
331  MaterialWrite(materialPtr);
332 }