ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4GDMLWriteParamvol.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4GDMLWriteParamvol.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 G4GDMLParamVol Implementation
29 //
30 // Original author: Zoltan Torzsok, November 2007
31 //
32 // --------------------------------------------------------------------
33 
34 #include "G4GDMLWriteParamvol.hh"
35 #include "G4GDMLWriteSolids.hh"
36 
37 #include "G4SystemOfUnits.hh"
38 #include "G4Box.hh"
39 #include "G4Trd.hh"
40 #include "G4Trap.hh"
41 #include "G4Tubs.hh"
42 #include "G4Cons.hh"
43 #include "G4Sphere.hh"
44 #include "G4Orb.hh"
45 #include "G4Torus.hh"
46 #include "G4Ellipsoid.hh"
47 #include "G4Para.hh"
48 #include "G4Hype.hh"
49 #include "G4Polycone.hh"
50 #include "G4Polyhedra.hh"
51 #include "G4LogicalVolume.hh"
52 #include "G4VPhysicalVolume.hh"
53 #include "G4PVParameterised.hh"
54 #include "G4VPVParameterisation.hh"
55 
58 {
59 }
60 
63 {
64 }
65 
67 Box_dimensionsWrite(xercesc::DOMElement* parametersElement,
68  const G4Box* const box)
69 {
70  xercesc::DOMElement* box_dimensionsElement = NewElement("box_dimensions");
71  box_dimensionsElement->
72  setAttributeNode(NewAttribute("x",2.0*box->GetXHalfLength()/mm));
73  box_dimensionsElement->
74  setAttributeNode(NewAttribute("y",2.0*box->GetYHalfLength()/mm));
75  box_dimensionsElement->
76  setAttributeNode(NewAttribute("z",2.0*box->GetZHalfLength()/mm));
77  box_dimensionsElement->
78  setAttributeNode(NewAttribute("lunit","mm"));
79  parametersElement->appendChild(box_dimensionsElement);
80 }
81 
83 Trd_dimensionsWrite(xercesc::DOMElement* parametersElement,
84  const G4Trd* const trd)
85 {
86  xercesc::DOMElement* trd_dimensionsElement = NewElement("trd_dimensions");
87  trd_dimensionsElement->
88  setAttributeNode(NewAttribute("x1",2.0*trd->GetXHalfLength1()/mm));
89  trd_dimensionsElement->
90  setAttributeNode(NewAttribute("x2",2.0*trd->GetXHalfLength2()/mm));
91  trd_dimensionsElement->
92  setAttributeNode(NewAttribute("y1",2.0*trd->GetYHalfLength1()/mm));
93  trd_dimensionsElement->
94  setAttributeNode(NewAttribute("y2",2.0*trd->GetYHalfLength2()/mm));
95  trd_dimensionsElement->
96  setAttributeNode(NewAttribute("z",2.0*trd->GetZHalfLength()/mm));
97  trd_dimensionsElement->
98  setAttributeNode(NewAttribute("lunit","mm"));
99  parametersElement->appendChild(trd_dimensionsElement);
100 }
101 
103 Trap_dimensionsWrite(xercesc::DOMElement* parametersElement,
104  const G4Trap* const trap)
105 {
106  const G4ThreeVector simaxis = trap->GetSymAxis();
107  const G4double phi = (simaxis.z() != 1.0)
108  ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
109  const G4double theta = std::acos(simaxis.z());
110  const G4double alpha1 = std::atan(trap->GetTanAlpha1());
111  const G4double alpha2 = std::atan(trap->GetTanAlpha2());
112 
113  xercesc::DOMElement* trap_dimensionsElement = NewElement("trap");
114  trap_dimensionsElement->
115  setAttributeNode(NewAttribute("z",2.0*trap->GetZHalfLength()/mm));
116  trap_dimensionsElement->
117  setAttributeNode(NewAttribute("theta",theta/degree));
118  trap_dimensionsElement->
119  setAttributeNode(NewAttribute("phi",phi/degree));
120  trap_dimensionsElement->
121  setAttributeNode(NewAttribute("y1",2.0*trap->GetYHalfLength1()/mm));
122  trap_dimensionsElement->
123  setAttributeNode(NewAttribute("x1",2.0*trap->GetXHalfLength1()/mm));
124  trap_dimensionsElement->
125  setAttributeNode(NewAttribute("x2",2.0*trap->GetXHalfLength2()/mm));
126  trap_dimensionsElement->
127  setAttributeNode(NewAttribute("alpha1",alpha1/degree));
128  trap_dimensionsElement->
129  setAttributeNode(NewAttribute("y2",2.0*trap->GetYHalfLength2()/mm));
130  trap_dimensionsElement->
131  setAttributeNode(NewAttribute("x3",2.0*trap->GetXHalfLength3()/mm));
132  trap_dimensionsElement->
133  setAttributeNode(NewAttribute("x4",2.0*trap->GetXHalfLength4()/mm));
134  trap_dimensionsElement->
135  setAttributeNode(NewAttribute("alpha2",alpha2/degree));
136  trap_dimensionsElement->
137  setAttributeNode(NewAttribute("aunit","deg"));
138  trap_dimensionsElement->
139  setAttributeNode(NewAttribute("lunit","mm"));
140  parametersElement->appendChild(trap_dimensionsElement);
141 }
142 
144 Tube_dimensionsWrite(xercesc::DOMElement* parametersElement,
145  const G4Tubs* const tube)
146 {
147  xercesc::DOMElement* tube_dimensionsElement = NewElement("tube_dimensions");
148  tube_dimensionsElement->
149  setAttributeNode(NewAttribute("InR",tube->GetInnerRadius()/mm));
150  tube_dimensionsElement->
151  setAttributeNode(NewAttribute("OutR",tube->GetOuterRadius()/mm));
152  tube_dimensionsElement->
153  setAttributeNode(NewAttribute("hz",2.0*tube->GetZHalfLength()/mm));
154  tube_dimensionsElement->
155  setAttributeNode(NewAttribute("StartPhi",tube->GetStartPhiAngle()/degree));
156  tube_dimensionsElement->
157  setAttributeNode(NewAttribute("DeltaPhi",tube->GetDeltaPhiAngle()/degree));
158  tube_dimensionsElement->
159  setAttributeNode(NewAttribute("aunit","deg"));
160  tube_dimensionsElement->
161  setAttributeNode(NewAttribute("lunit","mm"));
162  parametersElement->appendChild(tube_dimensionsElement);
163 }
164 
165 
167 Cone_dimensionsWrite(xercesc::DOMElement* parametersElement,
168  const G4Cons* const cone)
169 {
170  xercesc::DOMElement* cone_dimensionsElement = NewElement("cone_dimensions");
171  cone_dimensionsElement->
172  setAttributeNode(NewAttribute("rmin1",cone->GetInnerRadiusMinusZ()/mm));
173  cone_dimensionsElement->
174  setAttributeNode(NewAttribute("rmax1",cone->GetOuterRadiusMinusZ()/mm));
175  cone_dimensionsElement->
176  setAttributeNode(NewAttribute("rmin2",cone->GetInnerRadiusPlusZ()/mm));
177  cone_dimensionsElement->
178  setAttributeNode(NewAttribute("rmax2",cone->GetOuterRadiusPlusZ()/mm));
179  cone_dimensionsElement->
180  setAttributeNode(NewAttribute("z",2.0*cone->GetZHalfLength()/mm));
181  cone_dimensionsElement->
182  setAttributeNode(NewAttribute("startphi",cone->GetStartPhiAngle()/degree));
183  cone_dimensionsElement->
184  setAttributeNode(NewAttribute("deltaphi",cone->GetDeltaPhiAngle()/degree));
185  cone_dimensionsElement->
186  setAttributeNode(NewAttribute("aunit","deg"));
187  cone_dimensionsElement->
188  setAttributeNode(NewAttribute("lunit","mm"));
189  parametersElement->appendChild(cone_dimensionsElement);
190 }
191 
193 Sphere_dimensionsWrite(xercesc::DOMElement* parametersElement,
194  const G4Sphere* const sphere)
195 {
196  xercesc::DOMElement* sphere_dimensionsElement =
197  NewElement("sphere_dimensions");
198  sphere_dimensionsElement->setAttributeNode(NewAttribute("rmin",
199  sphere->GetInnerRadius()/mm));
200  sphere_dimensionsElement->setAttributeNode(NewAttribute("rmax",
201  sphere->GetOuterRadius()/mm));
202  sphere_dimensionsElement->setAttributeNode(NewAttribute("startphi",
203  sphere->GetStartPhiAngle()/degree));
204  sphere_dimensionsElement->setAttributeNode(NewAttribute("deltaphi",
205  sphere->GetDeltaPhiAngle()/degree));
206  sphere_dimensionsElement->setAttributeNode(NewAttribute("starttheta",
207  sphere->GetStartThetaAngle()/degree));
208  sphere_dimensionsElement->setAttributeNode(NewAttribute("deltatheta",
209  sphere->GetDeltaThetaAngle()/degree));
210  sphere_dimensionsElement->setAttributeNode(NewAttribute("aunit","deg"));
211  sphere_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
212  parametersElement->appendChild(sphere_dimensionsElement);
213 }
214 
216 Orb_dimensionsWrite(xercesc::DOMElement* parametersElement,
217  const G4Orb* const orb)
218 {
219  xercesc::DOMElement* orb_dimensionsElement = NewElement("orb_dimensions");
220  orb_dimensionsElement->setAttributeNode(NewAttribute("r",
221  orb->GetRadius()/mm));
222  orb_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
223  parametersElement->appendChild(orb_dimensionsElement);
224 }
225 
227 Torus_dimensionsWrite(xercesc::DOMElement* parametersElement,
228  const G4Torus* const torus)
229 {
230  xercesc::DOMElement* torus_dimensionsElement =
231  NewElement("torus_dimensions");
232  torus_dimensionsElement->
233  setAttributeNode(NewAttribute("rmin",torus->GetRmin()/mm));
234  torus_dimensionsElement->
235  setAttributeNode(NewAttribute("rmax",torus->GetRmax()/mm));
236  torus_dimensionsElement->
237  setAttributeNode(NewAttribute("rtor",torus->GetRtor()/mm));
238  torus_dimensionsElement->
239  setAttributeNode(NewAttribute("startphi",torus->GetSPhi()/degree));
240  torus_dimensionsElement->
241  setAttributeNode(NewAttribute("deltaphi",torus->GetDPhi()/degree));
242  torus_dimensionsElement->
243  setAttributeNode(NewAttribute("aunit","deg"));
244  torus_dimensionsElement->
245  setAttributeNode(NewAttribute("lunit","mm"));
246  parametersElement->appendChild(torus_dimensionsElement);
247 }
248 
250 Ellipsoid_dimensionsWrite(xercesc::DOMElement* parametersElement,
251  const G4Ellipsoid* const ellipsoid)
252 {
253  xercesc::DOMElement* ellipsoid_dimensionsElement =
254  NewElement("ellipsoid_dimensions");
255  ellipsoid_dimensionsElement->
256  setAttributeNode(NewAttribute("ax",ellipsoid->GetSemiAxisMax(0)/mm));
257  ellipsoid_dimensionsElement->
258  setAttributeNode(NewAttribute("by",ellipsoid->GetSemiAxisMax(1)/mm));
259  ellipsoid_dimensionsElement->
260  setAttributeNode(NewAttribute("cz",ellipsoid->GetSemiAxisMax(2)/mm));
261  ellipsoid_dimensionsElement->
262  setAttributeNode(NewAttribute("zcut1",ellipsoid->GetZBottomCut()/mm));
263  ellipsoid_dimensionsElement->
264  setAttributeNode(NewAttribute("zcut2",ellipsoid->GetZTopCut()/mm));
265  ellipsoid_dimensionsElement->
266  setAttributeNode(NewAttribute("lunit","mm"));
267  parametersElement->appendChild(ellipsoid_dimensionsElement);
268 }
269 
271 Para_dimensionsWrite(xercesc::DOMElement* parametersElement,
272  const G4Para* const para)
273 {
274  const G4ThreeVector simaxis = para->GetSymAxis();
275 
276  const G4double alpha = std::atan(para->GetTanAlpha());
277  const G4double theta = std::acos(simaxis.z());
278  const G4double phi = (simaxis.z() != 1.0)
279  ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
280 
281  xercesc::DOMElement* para_dimensionsElement = NewElement("para_dimensions");
282  para_dimensionsElement->
283  setAttributeNode(NewAttribute("x",2.0*para->GetXHalfLength()/mm));
284  para_dimensionsElement->
285  setAttributeNode(NewAttribute("y",2.0*para->GetYHalfLength()/mm));
286  para_dimensionsElement->
287  setAttributeNode(NewAttribute("z",2.0*para->GetZHalfLength()/mm));
288  para_dimensionsElement->
289  setAttributeNode(NewAttribute("alpha",alpha/degree));
290  para_dimensionsElement->
291  setAttributeNode(NewAttribute("theta",theta/degree));
292  para_dimensionsElement->
293  setAttributeNode(NewAttribute("phi",phi/degree));
294  para_dimensionsElement->
295  setAttributeNode(NewAttribute("aunit","deg"));
296  para_dimensionsElement->
297  setAttributeNode(NewAttribute("lunit","mm"));
298  parametersElement->appendChild(para_dimensionsElement);
299 }
300 
302 Hype_dimensionsWrite(xercesc::DOMElement* parametersElement,
303  const G4Hype* const hype)
304 {
305  xercesc::DOMElement* hype_dimensionsElement = NewElement("hype_dimensions");
306  hype_dimensionsElement->
307  setAttributeNode(NewAttribute("rmin",hype->GetInnerRadius()/mm));
308  hype_dimensionsElement->
309  setAttributeNode(NewAttribute("rmax",hype->GetOuterRadius()/mm));
310  hype_dimensionsElement->
311  setAttributeNode(NewAttribute("inst",hype->GetInnerStereo()/degree));
312  hype_dimensionsElement->
313  setAttributeNode(NewAttribute("outst",hype->GetOuterStereo()/degree));
314  hype_dimensionsElement->
315  setAttributeNode(NewAttribute("z",2.0*hype->GetZHalfLength()/mm));
316  hype_dimensionsElement->
317  setAttributeNode(NewAttribute("aunit","deg"));
318  hype_dimensionsElement->
319  setAttributeNode(NewAttribute("lunit","mm"));
320  parametersElement->appendChild(hype_dimensionsElement);
321 }
322 
324 Polycone_dimensionsWrite(xercesc::DOMElement* parametersElement,
325  const G4Polycone* const pcone)
326 {
327  xercesc::DOMElement* pcone_dimensionsElement
328  = NewElement("polycone_dimensions");
329 
330  pcone_dimensionsElement->setAttributeNode(NewAttribute("numRZ",
332  pcone_dimensionsElement->setAttributeNode(NewAttribute("startPhi",
334  pcone_dimensionsElement->setAttributeNode(NewAttribute("openPhi",
336  pcone_dimensionsElement->setAttributeNode(NewAttribute("aunit","deg"));
337  pcone_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
338 
339  parametersElement->appendChild(pcone_dimensionsElement);
340  const size_t num_zplanes = pcone->GetOriginalParameters()->Num_z_planes;
341  const G4double* z_array = pcone->GetOriginalParameters()->Z_values;
342  const G4double* rmin_array = pcone->GetOriginalParameters()->Rmin;
343  const G4double* rmax_array = pcone->GetOriginalParameters()->Rmax;
344 
345  for (size_t i=0; i<num_zplanes; i++)
346  {
347  ZplaneWrite(pcone_dimensionsElement,z_array[i],
348  rmin_array[i],rmax_array[i]);
349  }
350 }
351 
353 Polyhedra_dimensionsWrite(xercesc::DOMElement* parametersElement,
354  const G4Polyhedra* const polyhedra)
355 {
356  xercesc::DOMElement* polyhedra_dimensionsElement
357  = NewElement("polyhedra_dimensions");
358 
359  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("numRZ",
360  polyhedra->GetOriginalParameters()->Num_z_planes));
361  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("numSide",
362  polyhedra->GetOriginalParameters()->numSide));
363  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("startPhi",
364  polyhedra->GetOriginalParameters()->Start_angle/degree));
365  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("openPhi",
367  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("aunit","deg"));
368  polyhedra_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
369 
370  parametersElement->appendChild(polyhedra_dimensionsElement);
371  const size_t num_zplanes = polyhedra->GetOriginalParameters()->Num_z_planes;
372  const G4double* z_array = polyhedra->GetOriginalParameters()->Z_values;
373  const G4double* rmin_array = polyhedra->GetOriginalParameters()->Rmin;
374  const G4double* rmax_array = polyhedra->GetOriginalParameters()->Rmax;
375 
376  for (size_t i=0; i<num_zplanes; i++)
377  {
378  ZplaneWrite(polyhedra_dimensionsElement,z_array[i],
379  rmin_array[i],rmax_array[i]);
380  }
381 }
382 
384 ParametersWrite(xercesc::DOMElement* paramvolElement,
385  const G4VPhysicalVolume* const paramvol,const G4int& index)
386 {
387  paramvol->GetParameterisation()
388  ->ComputeTransformation(index, const_cast<G4VPhysicalVolume*>(paramvol));
389  G4ThreeVector Angles;
390  G4String name = GenerateName(paramvol->GetName(),paramvol);
391  std::stringstream os;
392  os.precision(15);
393  os << index;
394  G4String sncopie = os.str();
395 
396  xercesc::DOMElement* parametersElement = NewElement("parameters");
397  parametersElement->setAttributeNode(NewAttribute("number",index+1));
398 
399  PositionWrite(parametersElement, name+sncopie+"_pos",
400  paramvol->GetObjectTranslation());
401  Angles=GetAngles(paramvol->GetObjectRotationValue());
402  if (Angles.mag2()>DBL_EPSILON)
403  {
404  RotationWrite(parametersElement, name+sncopie+"_rot",
405  GetAngles(paramvol->GetObjectRotationValue()));
406  }
407  paramvolElement->appendChild(parametersElement);
408 
409  G4VSolid* solid = paramvol->GetLogicalVolume()->GetSolid();
410 
411  if (G4Box* box = dynamic_cast<G4Box*>(solid))
412  {
413  paramvol->GetParameterisation()->ComputeDimensions(*box,index,
414  const_cast<G4VPhysicalVolume*>(paramvol));
415  Box_dimensionsWrite(parametersElement,box);
416  } else
417  if (G4Trd* trd = dynamic_cast<G4Trd*>(solid))
418  {
419  paramvol->GetParameterisation()->ComputeDimensions(*trd,index,
420  const_cast<G4VPhysicalVolume*>(paramvol));
421  Trd_dimensionsWrite(parametersElement,trd);
422  } else
423  if (G4Trap* trap = dynamic_cast<G4Trap*>(solid))
424  {
425  paramvol->GetParameterisation()->ComputeDimensions(*trap,index,
426  const_cast<G4VPhysicalVolume*>(paramvol));
427  Trap_dimensionsWrite(parametersElement,trap);
428  } else
429  if (G4Tubs* tube = dynamic_cast<G4Tubs*>(solid))
430  {
431  paramvol->GetParameterisation()->ComputeDimensions(*tube,index,
432  const_cast<G4VPhysicalVolume*>(paramvol));
433  Tube_dimensionsWrite(parametersElement,tube);
434  } else
435  if (G4Cons* cone = dynamic_cast<G4Cons*>(solid))
436  {
437  paramvol->GetParameterisation()->ComputeDimensions(*cone,index,
438  const_cast<G4VPhysicalVolume*>(paramvol));
439  Cone_dimensionsWrite(parametersElement,cone);
440  } else
441  if (G4Sphere* sphere = dynamic_cast<G4Sphere*>(solid))
442  {
443  paramvol->GetParameterisation()->ComputeDimensions(*sphere,index,
444  const_cast<G4VPhysicalVolume*>(paramvol));
445  Sphere_dimensionsWrite(parametersElement,sphere);
446  } else
447  if (G4Orb* orb = dynamic_cast<G4Orb*>(solid))
448  {
449  paramvol->GetParameterisation()->ComputeDimensions(*orb,index,
450  const_cast<G4VPhysicalVolume*>(paramvol));
451  Orb_dimensionsWrite(parametersElement,orb);
452  } else
453  if (G4Torus* torus = dynamic_cast<G4Torus*>(solid))
454  {
455  paramvol->GetParameterisation()->ComputeDimensions(*torus,index,
456  const_cast<G4VPhysicalVolume*>(paramvol));
457  Torus_dimensionsWrite(parametersElement,torus);
458  } else
459  if (G4Ellipsoid* ellipsoid = dynamic_cast<G4Ellipsoid*>(solid))
460  {
461  paramvol->GetParameterisation()->ComputeDimensions(*ellipsoid,index,
462  const_cast<G4VPhysicalVolume*>(paramvol));
463  Ellipsoid_dimensionsWrite(parametersElement,ellipsoid);
464  } else
465  if (G4Para* para = dynamic_cast<G4Para*>(solid))
466  {
467  paramvol->GetParameterisation()->ComputeDimensions(*para,index,
468  const_cast<G4VPhysicalVolume*>(paramvol));
469  Para_dimensionsWrite(parametersElement,para);
470  } else
471  if (G4Hype* hype = dynamic_cast<G4Hype*>(solid))
472  {
473  paramvol->GetParameterisation()->ComputeDimensions(*hype,index,
474  const_cast<G4VPhysicalVolume*>(paramvol));
475  Hype_dimensionsWrite(parametersElement,hype);
476  }else
477  if (G4Polycone* pcone = dynamic_cast<G4Polycone*>(solid))
478  {
479  paramvol->GetParameterisation()->ComputeDimensions(*pcone,index,
480  const_cast<G4VPhysicalVolume*>(paramvol));
481  Polycone_dimensionsWrite(parametersElement,pcone);
482  }else
483  if (G4Polyhedra* polyhedra = dynamic_cast<G4Polyhedra*>(solid))
484  {
485  paramvol->GetParameterisation()->ComputeDimensions(*polyhedra,index,
486  const_cast<G4VPhysicalVolume*>(paramvol));
487  Polyhedra_dimensionsWrite(parametersElement,polyhedra);
488  }
489  else
490  {
491  G4String error_msg = "Solid '" + solid->GetName()
492  + "' cannot be used in parameterised volume!";
493  G4Exception("G4GDMLWriteParamvol::ParametersWrite()",
494  "InvalidSetup", FatalException, error_msg);
495  }
496 }
497 
499 ParamvolWrite(xercesc::DOMElement* volumeElement,
500  const G4VPhysicalVolume* const paramvol)
501 {
502  const G4String volumeref =
503  GenerateName(paramvol->GetLogicalVolume()->GetName(),
504  paramvol->GetLogicalVolume());
505  xercesc::DOMElement* paramvolElement = NewElement("paramvol");
506  paramvolElement->setAttributeNode(NewAttribute("ncopies",
507  paramvol->GetMultiplicity()));
508  xercesc::DOMElement* volumerefElement = NewElement("volumeref");
509  volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
510 
511  xercesc::DOMElement* algorithmElement =
512  NewElement("parameterised_position_size");
513  paramvolElement->appendChild(volumerefElement);
514  paramvolElement->appendChild(algorithmElement);
515  ParamvolAlgorithmWrite(algorithmElement,paramvol);
516  volumeElement->appendChild(paramvolElement);
517 }
518 
520 ParamvolAlgorithmWrite(xercesc::DOMElement* paramvolElement,
521  const G4VPhysicalVolume* const paramvol)
522 {
523  const G4String volumeref =
524  GenerateName(paramvol->GetLogicalVolume()->GetName(),
525  paramvol->GetLogicalVolume());
526 
527  const G4int parameterCount = paramvol->GetMultiplicity();
528 
529  for (G4int i=0; i<parameterCount; i++)
530  {
531  ParametersWrite(paramvolElement,paramvol,i);
532  }
533 }