G4GDMLWriteParamvol.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 //
00027 // $Id$
00028 //
00029 // class G4GDMLParamVol Implementation
00030 //
00031 // Original author: Zoltan Torzsok, November 2007
00032 //
00033 // --------------------------------------------------------------------
00034 
00035 #include "G4GDMLWriteParamvol.hh"
00036 
00037 #include "G4SystemOfUnits.hh"
00038 #include "G4Box.hh"
00039 #include "G4Trd.hh"
00040 #include "G4Trap.hh"
00041 #include "G4Tubs.hh"
00042 #include "G4Cons.hh"
00043 #include "G4Sphere.hh"
00044 #include "G4Orb.hh"
00045 #include "G4Torus.hh"
00046 #include "G4Para.hh"
00047 #include "G4Hype.hh"
00048 #include "G4LogicalVolume.hh"
00049 #include "G4VPhysicalVolume.hh"
00050 #include "G4PVParameterised.hh"
00051 #include "G4VPVParameterisation.hh"
00052 
00053 G4GDMLWriteParamvol::
00054 G4GDMLWriteParamvol() : G4GDMLWriteSetup()
00055 {
00056 }
00057 
00058 G4GDMLWriteParamvol::
00059 ~G4GDMLWriteParamvol()
00060 {
00061 }
00062 
00063 void G4GDMLWriteParamvol::
00064 Box_dimensionsWrite(xercesc::DOMElement* parametersElement,
00065                     const G4Box* const box)
00066 {
00067    xercesc::DOMElement* box_dimensionsElement = NewElement("box_dimensions");
00068    box_dimensionsElement->
00069      setAttributeNode(NewAttribute("x",2.0*box->GetXHalfLength()/mm));
00070    box_dimensionsElement->
00071      setAttributeNode(NewAttribute("y",2.0*box->GetYHalfLength()/mm));
00072    box_dimensionsElement->
00073      setAttributeNode(NewAttribute("z",2.0*box->GetZHalfLength()/mm));
00074    box_dimensionsElement->
00075      setAttributeNode(NewAttribute("lunit","mm"));
00076    parametersElement->appendChild(box_dimensionsElement);
00077 }
00078 
00079 void G4GDMLWriteParamvol::
00080 Trd_dimensionsWrite(xercesc::DOMElement* parametersElement,
00081                     const G4Trd* const trd)
00082 {
00083    xercesc::DOMElement* trd_dimensionsElement = NewElement("trd_dimensions");
00084    trd_dimensionsElement->
00085      setAttributeNode(NewAttribute("x1",2.0*trd->GetXHalfLength1()/mm));
00086    trd_dimensionsElement->
00087      setAttributeNode(NewAttribute("x2",2.0*trd->GetXHalfLength2()/mm));
00088    trd_dimensionsElement->
00089      setAttributeNode(NewAttribute("y1",2.0*trd->GetYHalfLength1()/mm));
00090    trd_dimensionsElement->
00091      setAttributeNode(NewAttribute("y2",2.0*trd->GetYHalfLength2()/mm));
00092    trd_dimensionsElement->
00093      setAttributeNode(NewAttribute("z",2.0*trd->GetZHalfLength()/mm));
00094    trd_dimensionsElement->
00095      setAttributeNode(NewAttribute("lunit","mm"));
00096    parametersElement->appendChild(trd_dimensionsElement);
00097 }
00098 
00099 void G4GDMLWriteParamvol::
00100 Trap_dimensionsWrite(xercesc::DOMElement* parametersElement,
00101                      const G4Trap* const trap)
00102 {
00103    const G4ThreeVector simaxis = trap->GetSymAxis();
00104    const G4double phi = (simaxis.z() != 1.0)
00105                       ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
00106    const G4double theta = std::acos(simaxis.z());
00107    const G4double alpha1 = std::atan(trap->GetTanAlpha1());
00108    const G4double alpha2 = std::atan(trap->GetTanAlpha2());
00109 
00110    xercesc::DOMElement* trap_dimensionsElement = NewElement("trap");
00111    trap_dimensionsElement->
00112      setAttributeNode(NewAttribute("z",2.0*trap->GetZHalfLength()/mm));
00113    trap_dimensionsElement->
00114      setAttributeNode(NewAttribute("theta",theta/degree));
00115    trap_dimensionsElement->
00116      setAttributeNode(NewAttribute("phi",phi/degree));
00117    trap_dimensionsElement->
00118      setAttributeNode(NewAttribute("y1",2.0*trap->GetYHalfLength1()/mm));
00119    trap_dimensionsElement->
00120      setAttributeNode(NewAttribute("x1",2.0*trap->GetXHalfLength1()/mm));
00121    trap_dimensionsElement->
00122      setAttributeNode(NewAttribute("x2",2.0*trap->GetXHalfLength2()/mm));
00123    trap_dimensionsElement->
00124      setAttributeNode(NewAttribute("alpha1",alpha1/degree));
00125    trap_dimensionsElement->
00126      setAttributeNode(NewAttribute("y2",2.0*trap->GetYHalfLength2()/mm));
00127    trap_dimensionsElement->
00128      setAttributeNode(NewAttribute("x3",2.0*trap->GetXHalfLength3()/mm));
00129    trap_dimensionsElement->
00130      setAttributeNode(NewAttribute("x4",2.0*trap->GetXHalfLength4()/mm));
00131    trap_dimensionsElement->
00132      setAttributeNode(NewAttribute("alpha2",alpha2/degree));
00133    trap_dimensionsElement->
00134      setAttributeNode(NewAttribute("aunit","deg"));
00135    trap_dimensionsElement->
00136      setAttributeNode(NewAttribute("lunit","mm"));
00137    parametersElement->appendChild(trap_dimensionsElement);
00138 }
00139 
00140 void G4GDMLWriteParamvol::
00141 Tube_dimensionsWrite(xercesc::DOMElement* parametersElement,
00142                      const G4Tubs* const tube)
00143 {
00144    xercesc::DOMElement* tube_dimensionsElement = NewElement("tube_dimensions");
00145    tube_dimensionsElement->
00146      setAttributeNode(NewAttribute("InR",tube->GetInnerRadius()/mm));
00147    tube_dimensionsElement->
00148      setAttributeNode(NewAttribute("OutR",tube->GetOuterRadius()/mm));
00149    tube_dimensionsElement->
00150      setAttributeNode(NewAttribute("hz",2.0*tube->GetZHalfLength()/mm));
00151    tube_dimensionsElement->
00152      setAttributeNode(NewAttribute("StartPhi",tube->GetStartPhiAngle()/degree));
00153    tube_dimensionsElement->
00154      setAttributeNode(NewAttribute("DeltaPhi",tube->GetDeltaPhiAngle()/degree));
00155    tube_dimensionsElement->
00156      setAttributeNode(NewAttribute("aunit","deg"));
00157    tube_dimensionsElement->
00158      setAttributeNode(NewAttribute("lunit","mm"));
00159    parametersElement->appendChild(tube_dimensionsElement);
00160 }
00161 
00162 
00163 void G4GDMLWriteParamvol::
00164 Cone_dimensionsWrite(xercesc::DOMElement* parametersElement,
00165                      const G4Cons* const cone)
00166 {
00167    xercesc::DOMElement* cone_dimensionsElement = NewElement("cone_dimensions");
00168    cone_dimensionsElement->
00169      setAttributeNode(NewAttribute("rmin1",cone->GetInnerRadiusMinusZ()/mm));
00170    cone_dimensionsElement->
00171      setAttributeNode(NewAttribute("rmax1",cone->GetOuterRadiusMinusZ()/mm));
00172    cone_dimensionsElement->
00173      setAttributeNode(NewAttribute("rmin2",cone->GetInnerRadiusPlusZ()/mm));
00174    cone_dimensionsElement->
00175      setAttributeNode(NewAttribute("rmax2",cone->GetOuterRadiusPlusZ()/mm));
00176    cone_dimensionsElement->
00177      setAttributeNode(NewAttribute("z",2.0*cone->GetZHalfLength()/mm));
00178    cone_dimensionsElement->
00179      setAttributeNode(NewAttribute("startphi",cone->GetStartPhiAngle()/degree));
00180    cone_dimensionsElement->
00181      setAttributeNode(NewAttribute("deltaphi",cone->GetDeltaPhiAngle()/degree));
00182    cone_dimensionsElement->
00183      setAttributeNode(NewAttribute("aunit","deg"));
00184    cone_dimensionsElement->
00185      setAttributeNode(NewAttribute("lunit","mm"));
00186    parametersElement->appendChild(cone_dimensionsElement);
00187 }
00188 
00189 void G4GDMLWriteParamvol::
00190 Sphere_dimensionsWrite(xercesc::DOMElement* parametersElement,
00191                        const G4Sphere* const sphere)
00192 {
00193    xercesc::DOMElement* sphere_dimensionsElement =
00194                         NewElement("sphere_dimensions");
00195    sphere_dimensionsElement->setAttributeNode(NewAttribute("rmin",
00196                              sphere->GetInsideRadius()/mm));
00197    sphere_dimensionsElement->setAttributeNode(NewAttribute("rmax",
00198                              sphere->GetOuterRadius()/mm));
00199    sphere_dimensionsElement->setAttributeNode(NewAttribute("startphi",
00200                              sphere->GetStartPhiAngle()/degree));
00201    sphere_dimensionsElement->setAttributeNode(NewAttribute("deltaphi",
00202                              sphere->GetDeltaPhiAngle()/degree));
00203    sphere_dimensionsElement->setAttributeNode(NewAttribute("starttheta",
00204                              sphere->GetStartThetaAngle()/degree));
00205    sphere_dimensionsElement->setAttributeNode(NewAttribute("deltatheta",
00206                              sphere->GetDeltaThetaAngle()/degree));
00207    sphere_dimensionsElement->setAttributeNode(NewAttribute("aunit","deg"));
00208    sphere_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
00209    parametersElement->appendChild(sphere_dimensionsElement);
00210 }
00211 
00212 void G4GDMLWriteParamvol::
00213 Orb_dimensionsWrite(xercesc::DOMElement* parametersElement,
00214                     const G4Orb* const orb)
00215 {
00216    xercesc::DOMElement* orb_dimensionsElement = NewElement("orb_dimensions");
00217    orb_dimensionsElement->setAttributeNode(NewAttribute("r",
00218                           orb->GetRadius()/mm));
00219    orb_dimensionsElement->setAttributeNode(NewAttribute("lunit","mm"));
00220    parametersElement->appendChild(orb_dimensionsElement);
00221 }
00222 
00223 void G4GDMLWriteParamvol::
00224 Torus_dimensionsWrite(xercesc::DOMElement* parametersElement,
00225                       const G4Torus* const torus)
00226 {
00227    xercesc::DOMElement* torus_dimensionsElement =
00228                         NewElement("torus_dimensions");
00229    torus_dimensionsElement->
00230      setAttributeNode(NewAttribute("rmin",torus->GetRmin()/mm));
00231    torus_dimensionsElement->
00232      setAttributeNode(NewAttribute("rmax",torus->GetRmax()/mm));
00233    torus_dimensionsElement->
00234      setAttributeNode(NewAttribute("rtor",torus->GetRtor()/mm));
00235    torus_dimensionsElement->
00236      setAttributeNode(NewAttribute("startphi",torus->GetSPhi()/degree));
00237    torus_dimensionsElement->
00238      setAttributeNode(NewAttribute("deltaphi",torus->GetDPhi()/degree));
00239    torus_dimensionsElement->
00240      setAttributeNode(NewAttribute("aunit","deg"));
00241    torus_dimensionsElement->
00242      setAttributeNode(NewAttribute("lunit","mm"));
00243    parametersElement->appendChild(torus_dimensionsElement);
00244 }
00245 
00246 void G4GDMLWriteParamvol::
00247 Para_dimensionsWrite(xercesc::DOMElement* parametersElement,
00248                      const G4Para* const para)
00249 {
00250    const G4ThreeVector simaxis = para->GetSymAxis();
00251    const G4double alpha = std::atan(para->GetTanAlpha());
00252    const G4double theta = std::acos(simaxis.z());
00253    const G4double phi = (simaxis.z() != 1.0)
00254                       ? (std::atan(simaxis.y()/simaxis.x())) : (0.0);
00255 
00256    xercesc::DOMElement* para_dimensionsElement = NewElement("para_dimensions");
00257    para_dimensionsElement->
00258      setAttributeNode(NewAttribute("x",2.0*para->GetXHalfLength()/mm));
00259    para_dimensionsElement->
00260      setAttributeNode(NewAttribute("y",2.0*para->GetYHalfLength()/mm));
00261    para_dimensionsElement->
00262      setAttributeNode(NewAttribute("z",2.0*para->GetZHalfLength()/mm));
00263    para_dimensionsElement->
00264      setAttributeNode(NewAttribute("alpha",alpha/degree));
00265    para_dimensionsElement->
00266      setAttributeNode(NewAttribute("theta",theta/degree));
00267    para_dimensionsElement->
00268      setAttributeNode(NewAttribute("phi",phi/degree));
00269    para_dimensionsElement->
00270      setAttributeNode(NewAttribute("aunit","deg"));
00271    para_dimensionsElement->
00272      setAttributeNode(NewAttribute("lunit","mm"));
00273    parametersElement->appendChild(para_dimensionsElement);
00274 }
00275 
00276 void G4GDMLWriteParamvol::
00277 Hype_dimensionsWrite(xercesc::DOMElement* parametersElement,
00278                      const G4Hype* const hype)
00279 {
00280    xercesc::DOMElement* hype_dimensionsElement = NewElement("hype_dimensions");
00281    hype_dimensionsElement->
00282      setAttributeNode(NewAttribute("rmin",hype->GetInnerRadius()/mm));
00283    hype_dimensionsElement->
00284      setAttributeNode(NewAttribute("rmax",hype->GetOuterRadius()/mm));
00285    hype_dimensionsElement->
00286      setAttributeNode(NewAttribute("inst",hype->GetInnerStereo()/degree));
00287    hype_dimensionsElement->
00288      setAttributeNode(NewAttribute("outst",hype->GetOuterStereo()/degree));
00289    hype_dimensionsElement->
00290      setAttributeNode(NewAttribute("z",2.0*hype->GetZHalfLength()/mm));
00291    hype_dimensionsElement->
00292      setAttributeNode(NewAttribute("aunit","deg"));
00293    hype_dimensionsElement->
00294      setAttributeNode(NewAttribute("lunit","mm"));
00295    parametersElement->appendChild(hype_dimensionsElement);
00296 }
00297 
00298 void G4GDMLWriteParamvol::
00299 ParametersWrite(xercesc::DOMElement* paramvolElement,
00300                 const G4VPhysicalVolume* const paramvol,const G4int& index)
00301 {
00302    paramvol->GetParameterisation()
00303      ->ComputeTransformation(index, const_cast<G4VPhysicalVolume*>(paramvol));
00304    G4ThreeVector Angles;
00305    G4String name = GenerateName(paramvol->GetName(),paramvol);
00306    std::stringstream os; 
00307    os.precision(15);
00308    os << index;     
00309    G4String sncopie = os.str(); 
00310 
00311    xercesc::DOMElement* parametersElement = NewElement("parameters");
00312    parametersElement->setAttributeNode(NewAttribute("number",index+1));
00313 
00314    PositionWrite(parametersElement, name+sncopie+"_pos",
00315                  paramvol->GetObjectTranslation());
00316    Angles=GetAngles(paramvol->GetObjectRotationValue());
00317    if (Angles.mag2()>DBL_EPSILON)
00318    {
00319      RotationWrite(parametersElement, name+sncopie+"_rot",
00320                    GetAngles(paramvol->GetObjectRotationValue()));
00321    }
00322    paramvolElement->appendChild(parametersElement);
00323 
00324    G4VSolid* solid = paramvol->GetLogicalVolume()->GetSolid();
00325 
00326    if (G4Box* box = dynamic_cast<G4Box*>(solid))
00327    {
00328       paramvol->GetParameterisation()->ComputeDimensions(*box,index,
00329                 const_cast<G4VPhysicalVolume*>(paramvol));
00330       Box_dimensionsWrite(parametersElement,box);
00331    } else
00332    if (G4Trd* trd = dynamic_cast<G4Trd*>(solid))
00333    {
00334       paramvol->GetParameterisation()->ComputeDimensions(*trd,index,
00335                 const_cast<G4VPhysicalVolume*>(paramvol));
00336       Trd_dimensionsWrite(parametersElement,trd);
00337    } else
00338    if (G4Trap* trap = dynamic_cast<G4Trap*>(solid))
00339    {
00340       paramvol->GetParameterisation()->ComputeDimensions(*trap,index,
00341                 const_cast<G4VPhysicalVolume*>(paramvol));
00342       Trap_dimensionsWrite(parametersElement,trap);
00343    } else
00344    if (G4Tubs* tube = dynamic_cast<G4Tubs*>(solid))
00345    {
00346       paramvol->GetParameterisation()->ComputeDimensions(*tube,index,
00347                 const_cast<G4VPhysicalVolume*>(paramvol));
00348       Tube_dimensionsWrite(parametersElement,tube);
00349    } else
00350    if (G4Cons* cone = dynamic_cast<G4Cons*>(solid))
00351    {
00352       paramvol->GetParameterisation()->ComputeDimensions(*cone,index,
00353                 const_cast<G4VPhysicalVolume*>(paramvol));
00354       Cone_dimensionsWrite(parametersElement,cone);
00355    } else
00356    if (G4Sphere* sphere = dynamic_cast<G4Sphere*>(solid))
00357    {
00358       paramvol->GetParameterisation()->ComputeDimensions(*sphere,index,
00359                 const_cast<G4VPhysicalVolume*>(paramvol));
00360       Sphere_dimensionsWrite(parametersElement,sphere);
00361    } else
00362    if (G4Orb* orb = dynamic_cast<G4Orb*>(solid))
00363    {
00364       paramvol->GetParameterisation()->ComputeDimensions(*orb,index,
00365                 const_cast<G4VPhysicalVolume*>(paramvol));
00366       Orb_dimensionsWrite(parametersElement,orb);
00367    } else
00368    if (G4Torus* torus = dynamic_cast<G4Torus*>(solid))
00369    {
00370       paramvol->GetParameterisation()->ComputeDimensions(*torus,index,
00371                 const_cast<G4VPhysicalVolume*>(paramvol));
00372       Torus_dimensionsWrite(parametersElement,torus);
00373    } else
00374    if (G4Para* para = dynamic_cast<G4Para*>(solid))
00375    {
00376       paramvol->GetParameterisation()->ComputeDimensions(*para,index,
00377                 const_cast<G4VPhysicalVolume*>(paramvol));
00378       Para_dimensionsWrite(parametersElement,para);
00379    } else
00380    if (G4Hype* hype = dynamic_cast<G4Hype*>(solid))
00381    {
00382       paramvol->GetParameterisation()->ComputeDimensions(*hype,index,
00383                 const_cast<G4VPhysicalVolume*>(paramvol));
00384       Hype_dimensionsWrite(parametersElement,hype);
00385    }
00386    else
00387    {
00388      G4String error_msg = "Solid '" + solid->GetName()
00389                         + "' cannot be used in parameterised volume!";
00390      G4Exception("G4GDMLWriteParamvol::ParametersWrite()",
00391                  "InvalidSetup", FatalException, error_msg);
00392    }
00393 }
00394 
00395 void G4GDMLWriteParamvol::
00396 ParamvolWrite(xercesc::DOMElement* volumeElement,
00397               const G4VPhysicalVolume* const paramvol)
00398 {
00399    const G4String volumeref =
00400                   GenerateName(paramvol->GetLogicalVolume()->GetName(),
00401                                paramvol->GetLogicalVolume());
00402    xercesc::DOMElement* paramvolElement = NewElement("paramvol");
00403    paramvolElement->setAttributeNode(NewAttribute("ncopies",
00404                                      paramvol->GetMultiplicity()));
00405    xercesc::DOMElement* volumerefElement = NewElement("volumeref");
00406    volumerefElement->setAttributeNode(NewAttribute("ref",volumeref));
00407 
00408    xercesc::DOMElement* algorithmElement =
00409      NewElement("parameterised_position_size");
00410    paramvolElement->appendChild(volumerefElement);
00411    paramvolElement->appendChild(algorithmElement);
00412    ParamvolAlgorithmWrite(algorithmElement,paramvol);
00413    volumeElement->appendChild(paramvolElement);
00414 }
00415 
00416 void G4GDMLWriteParamvol::
00417 ParamvolAlgorithmWrite(xercesc::DOMElement* paramvolElement,
00418                        const G4VPhysicalVolume* const paramvol)
00419 {
00420    const G4String volumeref =
00421                   GenerateName(paramvol->GetLogicalVolume()->GetName(),
00422                                paramvol->GetLogicalVolume());
00423   
00424    const G4int parameterCount = paramvol->GetMultiplicity();
00425 
00426    for (G4int i=0; i<parameterCount; i++)
00427    {
00428      ParametersWrite(paramvolElement,paramvol,i);
00429    }
00430 }

Generated on Mon May 27 17:48:20 2013 for Geant4 by  doxygen 1.4.7