Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4GDMLWriteMaterials.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4GDMLWriteMaterials.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 // $Id: PHG4GDMLWriteMaterials.cc 70764 2013-06-05 12:54:37Z gcosmo $
28 //
29 // class PHG4GDMLWriteMaterials Implementation
30 //
31 // Original author: Zoltan Torzsok, November 2007
32 //
33 // --------------------------------------------------------------------
34 
36 
37 #include <Geant4/G4Element.hh>
38 #include <Geant4/G4Isotope.hh>
39 #include <Geant4/G4Material.hh>
40 #include <Geant4/G4MaterialPropertiesTable.hh>
41 #include <Geant4/G4PhysicalConstants.hh>
42 #include <Geant4/G4PhysicsFreeVector.hh>
43 #include <Geant4/G4SystemOfUnits.hh>
44 
45 #include <sstream>
46 
49  , materialsElement(0)
50 {
51 }
52 
54 {
55 }
56 
58  AtomWrite(xercesc::DOMElement* element, const G4double& a)
59 {
60  xercesc::DOMElement* atomElement = NewElement("atom");
61  atomElement->setAttributeNode(NewAttribute("unit", "g/mole"));
62  atomElement->setAttributeNode(NewAttribute("value", a * mole / g));
63  element->appendChild(atomElement);
64 }
65 
67  DWrite(xercesc::DOMElement* element, const G4double& d)
68 {
69  xercesc::DOMElement* DElement = NewElement("D");
70  DElement->setAttributeNode(NewAttribute("unit", "g/cm3"));
71  DElement->setAttributeNode(NewAttribute("value", d * cm3 / g));
72  element->appendChild(DElement);
73 }
74 
76  PWrite(xercesc::DOMElement* element, const G4double& P)
77 {
78  xercesc::DOMElement* PElement = NewElement("P");
79  PElement->setAttributeNode(NewAttribute("unit", "pascal"));
80  PElement->setAttributeNode(NewAttribute("value", P / hep_pascal));
81  element->appendChild(PElement);
82 }
83 
85  TWrite(xercesc::DOMElement* element, const G4double& T)
86 {
87  xercesc::DOMElement* TElement = NewElement("T");
88  TElement->setAttributeNode(NewAttribute("unit", "K"));
89  TElement->setAttributeNode(NewAttribute("value", T / kelvin));
90  element->appendChild(TElement);
91 }
92 
94  MEEWrite(xercesc::DOMElement* element, const G4double& MEE)
95 {
96  xercesc::DOMElement* PElement = NewElement("MEE");
97  PElement->setAttributeNode(NewAttribute("unit", "eV"));
98  PElement->setAttributeNode(NewAttribute("value", MEE / electronvolt));
99  element->appendChild(PElement);
100 }
101 
103  IsotopeWrite(const G4Isotope* const isotopePtr)
104 {
105  const G4String name = GenerateName(isotopePtr->GetName(), isotopePtr);
106 
107  xercesc::DOMElement* isotopeElement = NewElement("isotope");
108  isotopeElement->setAttributeNode(NewAttribute("name", name));
109  isotopeElement->setAttributeNode(NewAttribute("N", isotopePtr->GetN()));
110  isotopeElement->setAttributeNode(NewAttribute("Z", isotopePtr->GetZ()));
111  materialsElement->appendChild(isotopeElement);
112  AtomWrite(isotopeElement, isotopePtr->GetA());
113 }
114 
115 void PHG4GDMLWriteMaterials::ElementWrite(const G4Element* const elementPtr)
116 {
117  const G4String name = GenerateName(elementPtr->GetName(), elementPtr);
118 
119  xercesc::DOMElement* elementElement = NewElement("element");
120  elementElement->setAttributeNode(NewAttribute("name", name));
121 
122  const size_t NumberOfIsotopes = elementPtr->GetNumberOfIsotopes();
123 
124  if (NumberOfIsotopes > 0)
125  {
126  const G4double* RelativeAbundanceVector =
127  elementPtr->GetRelativeAbundanceVector();
128  for (size_t i = 0; i < NumberOfIsotopes; i++)
129  {
130  G4String fractionref = GenerateName(elementPtr->GetIsotope(i)->GetName(),
131  elementPtr->GetIsotope(i));
132  xercesc::DOMElement* fractionElement = NewElement("fraction");
133  fractionElement->setAttributeNode(NewAttribute("n",
134  RelativeAbundanceVector[i]));
135  fractionElement->setAttributeNode(NewAttribute("ref", fractionref));
136  elementElement->appendChild(fractionElement);
137  AddIsotope(elementPtr->GetIsotope(i));
138  }
139  }
140  else
141  {
142  elementElement->setAttributeNode(NewAttribute("Z", elementPtr->GetZ()));
143  AtomWrite(elementElement, elementPtr->GetA());
144  }
145 
146  materialsElement->appendChild(elementElement);
147  // Append the element AFTER all the possible components are appended!
148 }
149 
150 void PHG4GDMLWriteMaterials::MaterialWrite(const G4Material* const materialPtr)
151 {
152  G4String state_str("undefined");
153  const G4State state = materialPtr->GetState();
154  if (state == kStateSolid)
155  {
156  state_str = "solid";
157  }
158  else if (state == kStateLiquid)
159  {
160  state_str = "liquid";
161  }
162  else if (state == kStateGas)
163  {
164  state_str = "gas";
165  }
166 
167  const G4String name = GenerateName(materialPtr->GetName(), materialPtr);
168 
169  xercesc::DOMElement* materialElement = NewElement("material");
170  materialElement->setAttributeNode(NewAttribute("name", name));
171  materialElement->setAttributeNode(NewAttribute("state", state_str));
172 
173  // Write any property attached to the material...
174  //
175  if (materialPtr->GetMaterialPropertiesTable())
176  {
177  PropertyWrite(materialElement, materialPtr);
178  }
179 
180  if (materialPtr->GetTemperature() != STP_Temperature)
181  {
182  TWrite(materialElement, materialPtr->GetTemperature());
183  }
184  if (materialPtr->GetPressure() != STP_Pressure)
185  {
186  PWrite(materialElement, materialPtr->GetPressure());
187  }
188 
189  // Write Ionisation potential (mean excitation energy)
190  MEEWrite(materialElement, materialPtr->GetIonisation()->GetMeanExcitationEnergy());
191 
192  DWrite(materialElement, materialPtr->GetDensity());
193 
194  const size_t NumberOfElements = materialPtr->GetNumberOfElements();
195 
196  if ((NumberOfElements > 1) || (materialPtr->GetElement(0) && materialPtr->GetElement(0)->GetNumberOfIsotopes() > 1))
197  {
198  const G4double* MassFractionVector = materialPtr->GetFractionVector();
199 
200  for (size_t i = 0; i < NumberOfElements; i++)
201  {
202  const G4String fractionref =
203  GenerateName(materialPtr->GetElement(i)->GetName(),
204  materialPtr->GetElement(i));
205  xercesc::DOMElement* fractionElement = NewElement("fraction");
206  fractionElement->setAttributeNode(NewAttribute("n",
207  MassFractionVector[i]));
208  fractionElement->setAttributeNode(NewAttribute("ref", fractionref));
209  materialElement->appendChild(fractionElement);
210  AddElement(materialPtr->GetElement(i));
211  }
212  }
213  else
214  {
215  materialElement->setAttributeNode(NewAttribute("Z", materialPtr->GetZ()));
216  AtomWrite(materialElement, materialPtr->GetA());
217  }
218 
219  // Append the material AFTER all the possible components are appended!
220  //
221  materialsElement->appendChild(materialElement);
222 }
223 
224 #if G4VERSION_NUMBER >= 1100
225 
226 // --------------------------------------------------------------------
227 void PHG4GDMLWriteMaterials::PropertyConstWrite(
228  const G4String& key, const G4double pval,
229  const G4MaterialPropertiesTable* ptable)
230 {
231  const G4String matrixref = GenerateName(key, ptable);
232  xercesc::DOMElement* matrixElement = NewElement("matrix");
233  matrixElement->setAttributeNode(NewAttribute("name", matrixref));
234  matrixElement->setAttributeNode(NewAttribute("coldim", "1"));
235  std::ostringstream pvalues;
236 
237  pvalues << pval;
238  matrixElement->setAttributeNode(NewAttribute("values", pvalues.str()));
239 
240  defineElement->appendChild(matrixElement);
241 }
242 
244  const G4String& key, const G4PhysicsFreeVector* const pvec)
245 {
246  for (std::size_t i = 0; i < propertyList.size(); ++i) // Check if property is
247  { // already in the list!
248  if (propertyList[i] == pvec)
249  {
250  return;
251  }
252  }
253  propertyList.push_back(pvec);
254 
255  const G4String matrixref = GenerateName(key, pvec);
256  xercesc::DOMElement* matrixElement = NewElement("matrix");
257  matrixElement->setAttributeNode(NewAttribute("name", matrixref));
258  matrixElement->setAttributeNode(NewAttribute("coldim", "2"));
259  std::ostringstream pvalues;
260  for (std::size_t i = 0; i < pvec->GetVectorLength(); ++i)
261  {
262  if (i != 0)
263  {
264  pvalues << " ";
265  }
266  pvalues << pvec->Energy(i) << " " << (*pvec)[i];
267  }
268  matrixElement->setAttributeNode(NewAttribute("values", pvalues.str()));
269 
270  defineElement->appendChild(matrixElement);
271 }
272 
273 void PHG4GDMLWriteMaterials::PropertyWrite(xercesc::DOMElement* matElement,
274  const G4Material* const mat)
275 {
276  xercesc::DOMElement* propElement;
277  G4MaterialPropertiesTable* ptable = mat->GetMaterialPropertiesTable();
278 
279  auto pvec = ptable->GetProperties();
280  auto cvec = ptable->GetConstProperties();
281 
282  for (size_t i = 0; i < pvec.size(); ++i)
283  {
284  if (pvec[i] != nullptr)
285  {
286  propElement = NewElement("property");
287  propElement->setAttributeNode(
288  NewAttribute("name", ptable->GetMaterialPropertyNames()[i]));
289  propElement->setAttributeNode(NewAttribute(
290  "ref", GenerateName(ptable->GetMaterialPropertyNames()[i],
291  pvec[i])));
292  PropertyVectorWrite(ptable->GetMaterialPropertyNames()[i],
293  pvec[i]);
294  matElement->appendChild(propElement);
295  }
296  }
297 
298  for (size_t i = 0; i < cvec.size(); ++i)
299  {
300  if (cvec[i].second == true)
301  {
302  propElement = NewElement("property");
303  propElement->setAttributeNode(NewAttribute(
304  "name", ptable->GetMaterialConstPropertyNames()[i]));
305  propElement->setAttributeNode(NewAttribute(
306  "ref", GenerateName(ptable->GetMaterialConstPropertyNames()[i],
307  ptable)));
308  PropertyConstWrite(ptable->GetMaterialConstPropertyNames()[i],
309  cvec[i].first, ptable);
310  matElement->appendChild(propElement);
311  }
312  }
313 }
314 
315 #else // #if G4VERSION_NUMBER >= 1100
317  const G4PhysicsOrderedFreeVector* const pvec)
318 {
319  const G4String matrixref = GenerateName(key, pvec);
320  xercesc::DOMElement* matrixElement = NewElement("matrix");
321  matrixElement->setAttributeNode(NewAttribute("name", matrixref));
322  matrixElement->setAttributeNode(NewAttribute("coldim", "2"));
323  std::ostringstream pvalues;
324  for (size_t i = 0; i < pvec->GetVectorLength(); i++)
325  {
326  if (i != 0)
327  {
328  pvalues << " ";
329  }
330  pvalues << pvec->Energy(i) << " " << (*pvec)[i];
331  }
332  matrixElement->setAttributeNode(NewAttribute("values", pvalues.str()));
333 
334  defineElement->appendChild(matrixElement);
335 }
336 
337 void PHG4GDMLWriteMaterials::PropertyWrite(xercesc::DOMElement* matElement,
338  const G4Material* const mat)
339 {
340  xercesc::DOMElement* propElement;
341  G4MaterialPropertiesTable* ptable = mat->GetMaterialPropertiesTable();
342 
343  const std::map<G4int, G4PhysicsOrderedFreeVector*,
344  std::less<G4int> >* pmap = ptable->GetPropertyMap();
345  const std::map<G4int, G4double,
346  std::less<G4int> >* cmap = ptable->GetConstPropertyMap();
347  std::map<G4int, G4PhysicsOrderedFreeVector*,
348  std::less<G4int> >::const_iterator mpos;
349  std::map<G4int, G4double,
350  std::less<G4int> >::const_iterator cpos;
351 
352  for (mpos = pmap->begin(); mpos != pmap->end(); ++mpos)
353  {
354  propElement = NewElement("property");
355  propElement->setAttributeNode(NewAttribute("name",
356  ptable->GetMaterialPropertyNames()[mpos->first]));
357  propElement->setAttributeNode(NewAttribute("ref",
358  GenerateName(ptable->GetMaterialPropertyNames()[mpos->first],
359  mpos->second)));
360  if (mpos->second)
361  {
362  PropertyVectorWrite(ptable->GetMaterialPropertyNames()[mpos->first],
363  mpos->second);
364  matElement->appendChild(propElement);
365  }
366  else
367  {
368  G4String warn_message = "Null pointer for material property -" + ptable->GetMaterialPropertyNames()[mpos->first] + "- of material -" + mat->GetName() + "- !";
369  G4Exception("G4GDMLWriteMaterials::PropertyWrite()", "NullPointer",
370  JustWarning, warn_message);
371  continue;
372  }
373  }
374 
375  for (cpos = cmap->begin(); cpos != cmap->end(); ++cpos)
376  {
377  propElement = NewElement("property");
378  propElement->setAttributeNode(NewAttribute("name",
379  ptable->GetMaterialConstPropertyNames()[cpos->first]));
380  propElement->setAttributeNode(NewAttribute("ref",
381  ptable->GetMaterialConstPropertyNames()[cpos->first]));
382  xercesc::DOMElement* constElement = NewElement("constant");
383  constElement->setAttributeNode(NewAttribute("name",
384  ptable->GetMaterialConstPropertyNames()[cpos->first]));
385  constElement->setAttributeNode(NewAttribute("value", cpos->second));
386  defineElement->appendChild(constElement);
387  matElement->appendChild(propElement);
388  }
389 }
390 #endif // #if G4VERSION_NUMBER >= 1100
391 
393 {
394  std::cout << "G4GDML: Writing materials..." << std::endl;
395 
396  materialsElement = NewElement("materials");
397  element->appendChild(materialsElement);
398 
399  isotopeList.clear();
400  elementList.clear();
401  materialList.clear();
402  propertyList.clear();
403 }
404 
405 void PHG4GDMLWriteMaterials::AddIsotope(const G4Isotope* const isotopePtr)
406 {
407  for (size_t i = 0; i < isotopeList.size(); i++) // Check if isotope is
408  { // already in the list!
409  if (isotopeList[i] == isotopePtr)
410  {
411  return;
412  }
413  }
414  isotopeList.push_back(isotopePtr);
415  IsotopeWrite(isotopePtr);
416 }
417 
418 void PHG4GDMLWriteMaterials::AddElement(const G4Element* const elementPtr)
419 {
420  for (size_t i = 0; i < elementList.size(); i++) // Check if element is
421  { // already in the list!
422  if (elementList[i] == elementPtr)
423  {
424  return;
425  }
426  }
427  elementList.push_back(elementPtr);
428  ElementWrite(elementPtr);
429 }
430 
431 void PHG4GDMLWriteMaterials::AddMaterial(const G4Material* const materialPtr)
432 {
433  for (size_t i = 0; i < materialList.size(); i++) // Check if material is
434  { // already in the list!
435  if (materialList[i] == materialPtr)
436  {
437  return;
438  }
439  }
440  materialList.push_back(materialPtr);
441  MaterialWrite(materialPtr);
442 }