Newer
Older
Markus Frank
committed
//==========================================================================
Markus Frank
committed
//--------------------------------------------------------------------------
// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
Markus Frank
committed
// All rights reserved.
Markus Frank
committed
// For the licensing terms see $DD4hepINSTALL/LICENSE.
// For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
Markus Frank
committed
// Author : M.Frank
//
//==========================================================================
Markus Frank
committed
// Framework include files
#include "DD4hep/Plugins.h"
#include "DD4hep/Shapes.h"
Markus Frank
committed
#include "DD4hep/Printout.h"
#include "DD4hep/DD4hepUnits.h"
#include "DD4hep/PropertyTable.h"
#include "DD4hep/detail/ShapesInterna.h"
#include "DD4hep/detail/ObjectsInterna.h"
#include "DD4hep/detail/DetectorInterna.h"
#include "DDG4/Geant4Field.h"
#include "DDG4/Geant4Converter.h"
#include "Geant4ShapeConverter.h"
// ROOT includes
#include "TROOT.h"
#include "TColor.h"
Markus Frank
committed
#include "TGeoNode.h"
#include "TGeoMatrix.h"
#include "TGeoBoolNode.h"
Markus Frank
committed
#include "TGeoParaboloid.h"
#include "TGeoManager.h"
#include "TClass.h"
#include "TMath.h"
#include "G4VisAttributes.hh"
#include "G4ProductionCuts.hh"
#include "G4VUserRegionInformation.hh"
#include "G4Ellipsoid.hh"
#include "G4UnionSolid.hh"
#include "G4ReflectedSolid.hh"
#include "G4EllipticalTube.hh"
#include "G4SubtractionSolid.hh"
#include "G4IntersectionSolid.hh"
#include "G4Element.hh"
#include "G4Region.hh"
#include "G4UserLimits.hh"
#include "G4LogicalVolume.hh"
#include "G4Material.hh"
#include "G4Element.hh"
#include "G4Isotope.hh"
#include "G4Transform3D.hh"
#include "G4ThreeVector.hh"
#include "G4PVPlacement.hh"
#include "G4ElectroMagneticField.hh"
#include "G4FieldManager.hh"
#include "G4ReflectionFactory.hh"
#include "G4OpticalSurface.hh"
#include "G4LogicalSkinSurface.hh"
#include "G4LogicalBorderSurface.hh"
#include "G4MaterialPropertiesTable.hh"
#include "G4MaterialPropertiesIndex.hh"
#if G4VERSION_NUMBER >= 1030
#include "G4ScaledSolid.hh"
#endif
#include "CLHEP/Units/SystemOfUnits.h"
Markus Frank
committed
// C/C++ include files
Markus Frank
committed
#include <iostream>
#include <iomanip>
#include <sstream>
namespace units = dd4hep;
using namespace dd4hep::detail;
using namespace dd4hep::sim;
using namespace dd4hep;
#include "DDG4/Geant4AssemblyVolume.h"
#include "DD4hep/DetectorTools.h"
#include "G4RotationMatrix.hh"
#include "G4AffineTransform.hh"
#include "G4LogicalVolume.hh"
#include "G4VPhysicalVolume.hh"
#include "G4ReflectionFactory.hh"
static const double CM_2_MM = (CLHEP::centimeter/dd4hep::centimeter);
void Geant4AssemblyVolume::imprint(Geant4GeometryInfo& info,
const TGeoNode* parent,
Chain chain,
Markus Frank
committed
Geant4AssemblyVolume* pAssembly,
G4LogicalVolume* pMotherLV,
G4Transform3D& transformation,
G4int copyNumBase,
G4bool surfCheck)
TGeoVolume* vol = parent->GetVolume();
unsigned int numberOfDaughters = (copyNumBase == 0) ? pMotherLV->GetNoDaughters() : copyNumBase;
++level;
// We start from the first available index
numberOfDaughters++;
ImprintsCountPlus();
Markus Frank
committed
vector<G4AssemblyTriplet> triplets = pAssembly->fTriplets;
//cout << " Assembly:" << detail::tools::placementPath(chain) << endl;
for( unsigned int i = 0; i < triplets.size(); i++ ) {
const TGeoNode* node = pAssembly->m_entries[i];
Chain new_chain = chain;
new_chain.emplace_back(node);
//cout << " Assembly: Entry: " << detail::tools::placementPath(new_chain) << endl;
G4Transform3D Ta( *(triplets[i].GetRotation()),
triplets[i].GetTranslation() );
if ( triplets[i].IsReflection() ) {
Ta = Ta * G4ReflectZ3D();
}
G4Transform3D Tfinal = transformation * Ta;
if ( triplets[i].GetVolume() ) {
Markus Frank
committed
// Generate the unique name for the next PV instance
// The name has format:
//
// av_WWW_impr_XXX_YYY_ZZZ
// where the fields mean:
// WWW - assembly volume instance number
// XXX - assembly volume imprint number
// YYY - the name of a log. volume we want to make a placement of
// ZZZ - the log. volume index inside the assembly volume
//
stringstream pvName;
Markus Frank
committed
pvName << "av_"
<< GetAssemblyID()
<< "_impr_"
<< GetImprintsCount()
<< "_"
<< triplets[i].GetVolume()->GetName().c_str()
<< "_pv_"
<< i
<< ends;
Markus Frank
committed
// Generate a new physical volume instance inside a mother
// (as we allow 3D transformation use G4ReflectionFactory to
// take into account eventual reflection)
//
#if 0
printout(INFO,"Geant4Converter","+++ Place %svolume %s in assembly.",
triplets[i].IsReflection() ? "REFLECTED " : "",
detail::tools::placementPath(new_chain).c_str());
#endif
Markus Frank
committed
G4PhysicalVolumesPair pvPlaced
= G4ReflectionFactory::Instance()->Place( Tfinal,
pvName.str().c_str(),
triplets[i].GetVolume(),
pMotherLV,
false,
numberOfDaughters + i,
surfCheck );
// Register the physical volume created by us so we can delete it later
//
//fPVStore.emplace_back( pvPlaced.first );
info.g4VolumeImprints[vol].emplace_back(new_chain,pvPlaced.first);
#if 0
Markus Frank
committed
cout << " Assembly:Parent:" << parent->GetName() << " " << node->GetName()
<< " " << (void*)node << " G4:" << pvName.str() << " Daughter:"
<< detail::tools::placementPath(new_chain) << endl;
Markus Frank
committed
cout << endl;
#endif
Markus Frank
committed
if ( pvPlaced.second ) {
G4Exception("G4AssemblyVolume::MakeImprint(..)", "GeomVol0003", FatalException,
"Fancy construct popping new mother from the stack!");
//fPVStore.emplace_back( pvPlaced.second );
Markus Frank
committed
}
else if ( triplets[i].GetAssembly() ) {
// Place volumes in this assembly with composed transformation
imprint(info, parent, new_chain, (Geant4AssemblyVolume*)triplets[i].GetAssembly(),
pMotherLV, Tfinal, i*100+copyNumBase, surfCheck );
else {
--level;
G4Exception("G4AssemblyVolume::MakeImprint(..)",
"GeomVol0003", FatalException,
"Triplet has no volume and no assembly");
Markus Frank
committed
}
//cout << "Imprinted assembly level:" << level << " in mother:" << pMotherLV->GetName() << endl;
Markus Frank
committed
}
namespace {
static string indent = "";
static Double_t s_identity_rot[] = { 1., 0., 0., 0., 1., 0., 0., 0., 1. };
struct MyTransform3D : public G4Transform3D {
MyTransform3D(double XX, double XY, double XZ, double DX, double YX, double YY, double YZ, double DY, double ZX, double ZY,
Markus Frank
committed
double ZZ, double DZ)
: G4Transform3D(XX, XY, XZ, DX, YX, YY, YZ, DY, ZX, ZY, ZZ, DZ) {
Markus Frank
committed
MyTransform3D(const double* t, const double* r = s_identity_rot)
: G4Transform3D(r[0],r[1],r[2],t[0]*CM_2_MM,r[3],r[4],r[5],t[1]*CM_2_MM,r[6],r[7],r[8],t[2]*CM_2_MM) {
}
MyTransform3D(Transform3D&& copy) : Transform3D(copy) {}
bool is_left_handed(const TGeoMatrix* m) {
const Double_t* r = m->GetRotationMatrix();
if ( r ) {
Double_t det =
r[0]*r[4]*r[8] + r[3]*r[7]*r[2] + r[6]*r[1]*r[5] -
r[2]*r[4]*r[6] - r[5]*r[7]*r[0] - r[8]*r[1]*r[3];
return det < 0e0;
}
return false;
}
class G4UserRegionInformation : public G4VUserRegionInformation {
public:
Region region;
double threshold;
bool storeSecondaries;
G4UserRegionInformation()
Markus Frank
committed
: threshold(0.0), storeSecondaries(false) {
}
virtual ~G4UserRegionInformation() {
}
if (region.isValid())
printout(DEBUG, "Region", "Name:%s", region.name());
pair<double,double> g4PropertyConversion(int index) {
switch(index) {
case kRINDEX: return make_pair(CLHEP::keV/units::keV, 1.0);
case kREFLECTIVITY: return make_pair(CLHEP::keV/units::keV, 1.0);
case kREALRINDEX: return make_pair(CLHEP::keV/units::keV, 1.0);
case kIMAGINARYRINDEX: return make_pair(CLHEP::keV/units::keV, 1.0);
case kEFFICIENCY: return make_pair(CLHEP::keV/units::keV, 1.0);
case kTRANSMITTANCE: return make_pair(CLHEP::keV/units::keV, 1.0);
case kSPECULARLOBECONSTANT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kSPECULARSPIKECONSTANT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kBACKSCATTERCONSTANT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kGROUPVEL: return make_pair(CLHEP::keV/units::keV, (CLHEP::m/CLHEP::s)/(units::m/units::s)); // meter/second
case kMIEHG: return make_pair(CLHEP::keV/units::keV, CLHEP::m/units::m);
case kRAYLEIGH: return make_pair(CLHEP::keV/units::keV, CLHEP::m/units::m); // ??? says its a length
case kWLSCOMPONENT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kWLSABSLENGTH: return make_pair(CLHEP::keV/units::keV, CLHEP::m/units::m);
case kABSLENGTH: return make_pair(CLHEP::keV/units::keV, CLHEP::m/units::m);
case kFASTCOMPONENT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kSLOWCOMPONENT: return make_pair(CLHEP::keV/units::keV, 1.0);
case kPROTONSCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV); // Yields: 1/energy
case kDEUTERONSCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV);
case kTRITONSCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV);
case kALPHASCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV);
case kIONSCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV);
case kELECTRONSCINTILLATIONYIELD: return make_pair(CLHEP::keV/units::keV, units::keV/CLHEP::keV);
default:
break;
}
printout(FATAL,"Geant4Converter", "+++ Cannot convert material property with index: %d", index);
#else
printout(FATAL,"Geant4Converter", "+++ Cannot convert material property with index: %d [Need Geant4 > 10.03]", index);
#endif
return make_pair(0e0,0e0);
}
double g4ConstPropertyConversion(int index) {
case kSURFACEROUGHNESS: return CLHEP::m/units::m; // Length
case kISOTHERMAL_COMPRESSIBILITY: return (CLHEP::m3/CLHEP::keV)/(units::m3/CLHEP::keV); // Volume/Energy
case kRS_SCALE_FACTOR: return 1.0; // ??
case kWLSMEANNUMBERPHOTONS: return 1.0; // ??
case kWLSTIMECONSTANT: return CLHEP::second/units::second; // Time
case kMIEHG_FORWARD: return 1.0;
case kMIEHG_BACKWARD: return 1.0;
case kMIEHG_FORWARD_RATIO: return 1.0;
case kSCINTILLATIONYIELD: return units::keV/CLHEP::keV; // Energy
case kRESOLUTIONSCALE: return 1.0;
case kFASTTIMECONSTANT: return CLHEP::second/units::second; // Time
case kFASTSCINTILLATIONRISETIME: return CLHEP::second/units::second; // Time
case kSLOWTIMECONSTANT: return CLHEP::second/units::second; // Time
case kSLOWSCINTILLATIONRISETIME: return CLHEP::second/units::second; // Time
case kFERMIPOT: return CLHEP::keV/units::keV; // Energy
case kDIFFUSION: return 1.0;
case kSPINFLIP: return 1.0;
case kLOSS: return 1.0; // ??
case kLOSSCS: return CLHEP::barn/units::barn; // ??
case kABSCS: return CLHEP::barn/units::barn; // ??
case kSCATCS: return CLHEP::barn/units::barn; // ??
case kMR_NBTHETA: return 1.0;
case kMR_NBE: return 1.0;
case kMR_RRMS: return 1.0; // ??
case kMR_CORRLEN: return CLHEP::m/units::m; // Length
case kMR_THETAMIN: return 1.0;
case kMR_THETAMAX: return 1.0;
case kMR_EMIN: return CLHEP::keV/units::keV; // Energy
case kMR_EMAX: return CLHEP::keV/units::keV; // Energy
case kMR_ANGNOTHETA: return 1.0;
case kMR_ANGNOPHI: return 1.0;
case kMR_ANGCUT: return 1.0;
default:
break;
}
printout(FATAL,"Geant4Converter", "+++ Cannot convert CONST material property with index: %d", index);
#else
printout(FATAL,"Geant4Converter", "+++ Cannot convert material property with index: %d [Need Geant4 > 10.03]", index);
#endif
Geant4Converter::Geant4Converter(const Detector& description_ref)
: Geant4Mapping(description_ref), checkOverlaps(true) {
this->Geant4Mapping::init();
Markus Frank
committed
m_propagateRegions = true;
Markus Frank
committed
outputLevel = PrintLevel(printLevel() - 1);
Geant4Converter::Geant4Converter(const Detector& description_ref, PrintLevel level)
: Geant4Mapping(description_ref), checkOverlaps(true) {
Markus Frank
committed
m_propagateRegions = true;
Markus Frank
committed
outputLevel = level;
Markus Frank
committed
/// Standard destructor
Geant4Converter::~Geant4Converter() {
Markus Frank
committed
/// Handle the conversion of isotopes
void* Geant4Converter::handleIsotope(const string& /* name */, const TGeoIsotope* iso) const {
G4Isotope* g4i = data().g4Isotopes[iso];
if (!g4i) {
double a_conv = (CLHEP::g / CLHEP::mole);
g4i = new G4Isotope(iso->GetName(), iso->GetZ(), iso->GetN(), iso->GetA()*a_conv);
printout(debugElements ? ALWAYS : outputLevel,
"Geant4Converter", "++ Created G4 Isotope %s from data: Z=%d N=%d A=%.3f [g/mole]",
iso->GetName(), iso->GetZ(), iso->GetN(), iso->GetA());
data().g4Isotopes[iso] = g4i;
}
return g4i;
}
/// Handle the conversion of elements
void* Geant4Converter::handleElement(const string& name, const Atom element) const {
G4Element* g4e = data().g4Elements[element];
PrintLevel lvl = debugElements ? ALWAYS : outputLevel;
if (element->GetNisotopes() > 0) {
g4e = new G4Element(name, element->GetTitle(), element->GetNisotopes());
for (int i = 0, n = element->GetNisotopes(); i < n; ++i) {
TGeoIsotope* iso = element->GetIsotope(i);
G4Isotope* g4iso = (G4Isotope*)handleIsotope(iso->GetName(), iso);
g4e->AddIsotope(g4iso, element->GetRelativeAbundance(i));
else {
// This adds in Geant4 the natural isotopes, which we normally do not want. We want to steer it outselves.
double a_conv = (CLHEP::g / CLHEP::mole);
g4e = new G4Element(element->GetTitle(), name, element->Z(), element->A()*a_conv);
printout(lvl, "Geant4Converter", "++ Created G4 Isotope %s from data: Z=%d N=%d A=%.3f [g/mole]",
element->GetName(), element->Z(), element->N(), element->A());
}
stringstream str;
str << (*g4e) << endl;
printout(lvl, "Geant4Converter", "++ Created G4 element %s", str.str().c_str());
data().g4Elements[element] = g4e;
}
return g4e;
}
/// Dump material in GDML format to output stream
void* Geant4Converter::handleMaterial(const string& name, Material medium) const {
Geant4GeometryInfo& info = data();
G4Material* mat = info.g4Materials[medium];
if ( !mat ) {
Markus Frank
committed
PrintLevel lvl = debugMaterials ? ALWAYS : outputLevel;
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
TGeoMaterial* material = medium->GetMaterial();
G4State state = kStateUndefined;
double density = material->GetDensity() * (CLHEP::gram / CLHEP::cm3);
if (density < 1e-25)
density = 1e-25;
switch (material->GetState()) {
case TGeoMaterial::kMatStateSolid:
state = kStateSolid;
break;
case TGeoMaterial::kMatStateLiquid:
state = kStateLiquid;
break;
case TGeoMaterial::kMatStateGas:
state = kStateGas;
break;
default:
case TGeoMaterial::kMatStateUndefined:
state = kStateUndefined;
break;
}
if (material->IsMixture()) {
double A_total = 0.0;
double W_total = 0.0;
TGeoMixture* mix = (TGeoMixture*) material;
int nElements = mix->GetNelements();
mat = new G4Material(name, density, nElements, state,
material->GetTemperature(), material->GetPressure());
for (int i = 0; i < nElements; ++i) {
A_total += (mix->GetAmixt())[i];
W_total += (mix->GetWmixt())[i];
for (int i = 0; i < nElements; ++i) {
TGeoElement* e = mix->GetElement(i);
G4Element* g4e = (G4Element*) handleElement(e->GetName(), Atom(e));
if (!g4e) {
printout(ERROR, "Material",
"Missing component %s for material %s. A=%f W=%f",
e->GetName(), mix->GetName(), A_total, W_total);
Markus Frank
committed
}
//mat->AddElement(g4e, (mix->GetAmixt())[i] / A_total);
mat->AddElement(g4e, (mix->GetWmixt())[i] / W_total);
}
else {
double z = material->GetZ(), a = material->GetA();
if ( z < 1.0000001 ) z = 1.0;
if ( a < 0.5000001 ) a = 1.0;
mat = new G4Material(name, z, a, density, state,
material->GetTemperature(), material->GetPressure());
}
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
/// Attach the material properties if any
G4MaterialPropertiesTable* tab = 0;
TListIter propIt(&material->GetProperties());
for(TObject* obj=propIt.Next(); obj; obj = propIt.Next()) {
TNamed* named = (TNamed*)obj;
TGDMLMatrix* matrix = info.manager->GetGDMLMatrix(named->GetTitle());
Geant4GeometryInfo::PropertyVector* v =
(Geant4GeometryInfo::PropertyVector*)handleMaterialProperties(matrix);
if ( 0 == v ) {
except("Geant4Converter", "++ FAILED to create G4 material %s [Cannot convert property:%s]",
material->GetName(), named->GetName());
if ( 0 == tab ) {
tab = new G4MaterialPropertiesTable();
mat->SetMaterialPropertiesTable(tab);
int idx = tab->GetPropertyIndex(named->GetName(), false);
if ( idx < 0 ) {
printout(ERROR, "Geant4Converter", "++ UNKNOWN Geant4 CONST Property: %-20s [IGNORED]", named->GetName());
continue;
}
// We need to convert the property from TGeo units to Geant4 units
auto conv = g4PropertyConversion(idx);
vector<double> bins(v->bins), vals(v->values);
for(size_t i=0, count=bins.size(); i<count; ++i)
bins[i] *= conv.first, vals[i] *= conv.second;
G4MaterialPropertyVector* vec = new G4MaterialPropertyVector(&bins[0], &vals[0], bins.size());
tab->AddProperty(named->GetName(), vec);
printout(lvl, "Geant4Converter", "++ Property: %-20s [%ld x %ld] -> %s ",
named->GetName(), matrix->GetRows(), matrix->GetCols(), named->GetTitle());
for(size_t i=0, count=v->bins.size(); i<count; ++i)
printout(lvl,named->GetName()," Geant4: %8.3g [MeV] TGeo: %8.3g [GeV] Conversion: %8.3g",
bins[i], v->bins[i], conv.first);
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
/// Attach the material properties if any
TListIter cpropIt(&material->GetConstProperties());
for(TObject* obj=cpropIt.Next(); obj; obj = cpropIt.Next()) {
Bool_t err = kFALSE;
TNamed* named = (TNamed*)obj;
double v = info.manager->GetProperty(named->GetTitle(),&err);
if ( err != kFALSE ) {
except("Geant4Converter",
"++ FAILED to create G4 material %s "
"[Cannot convert const property: %s]",
material->GetName(), named->GetName());
}
if ( 0 == tab ) {
tab = new G4MaterialPropertiesTable();
mat->SetMaterialPropertiesTable(tab);
}
int idx = tab->GetConstPropertyIndex(named->GetName(), false);
if ( idx < 0 ) {
printout(ERROR, "Geant4Converter", "++ UNKNOWN Geant4 CONST Property: %-20s [IGNORED]", named->GetName());
continue;
}
// We need to convert the property from TGeo units to Geant4 units
double conv = g4ConstPropertyConversion(idx);
printout(lvl, "Geant4Converter", "++ CONST Property: %-20s %g ", named->GetName(), v);
tab->AddConstProperty(named->GetName(), v * conv);
}
#endif
auto* ionization = mat->GetIonisation();
stringstream str;
str << (*mat) << endl;
if ( ionization )
str << " log(MEE): " << ionization->GetLogMeanExcEnergy();
else
str << " log(MEE): UNKNOWN";
printout(lvl, "Geant4Converter", "++ Created G4 material %s", str.str().c_str());
info.g4Materials[medium] = mat;
void* Geant4Converter::handleSolid(const string& name, const TGeoShape* shape) const {
Markus Frank
committed
G4VSolid* solid = 0;
if ( shape ) {
if ( 0 != (solid = data().g4Solids[shape]) ) {
Markus Frank
committed
return solid;
}
TClass* isa = shape->IsA();
PrintLevel lvl = debugShapes ? ALWAYS : outputLevel;
if (isa == TGeoShapeAssembly::Class()) {
// Assemblies have no corresponding 'shape' in Geant4. Ignore the shape translation.
// It does not harm, since this 'shape' is never accessed afterwards.
data().g4Solids[shape] = solid = convertShape<TGeoShapeAssembly>(shape);
return solid;
Markus Frank
committed
}
else if (isa == TGeoBBox::Class())
solid = convertShape<TGeoBBox>(shape);
else if (isa == TGeoTube::Class())
solid = convertShape<TGeoTube>(shape);
else if (isa == TGeoTubeSeg::Class())
solid = convertShape<TGeoTubeSeg>(shape);
else if (isa == TGeoCtub::Class())
solid = convertShape<TGeoCtub>(shape);
else if (isa == TGeoEltu::Class())
solid = convertShape<TGeoEltu>(shape);
else if (isa == TwistedTubeObject::Class())
solid = convertShape<TwistedTubeObject>(shape);
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
else if (isa == TGeoTrd1::Class())
solid = convertShape<TGeoTrd1>(shape);
else if (isa == TGeoTrd2::Class())
solid = convertShape<TGeoTrd2>(shape);
else if (isa == TGeoHype::Class())
solid = convertShape<TGeoHype>(shape);
else if (isa == TGeoXtru::Class())
solid = convertShape<TGeoXtru>(shape);
else if (isa == TGeoPgon::Class())
solid = convertShape<TGeoPgon>(shape);
else if (isa == TGeoPcon::Class())
solid = convertShape<TGeoPcon>(shape);
else if (isa == TGeoCone::Class())
solid = convertShape<TGeoCone>(shape);
else if (isa == TGeoConeSeg::Class())
solid = convertShape<TGeoConeSeg>(shape);
else if (isa == TGeoParaboloid::Class())
solid = convertShape<TGeoParaboloid>(shape);
else if (isa == TGeoSphere::Class())
solid = convertShape<TGeoSphere>(shape);
else if (isa == TGeoTorus::Class())
solid = convertShape<TGeoTorus>(shape);
else if (isa == TGeoTrap::Class())
solid = convertShape<TGeoTrap>(shape);
else if (isa == TGeoArb8::Class())
solid = convertShape<TGeoArb8>(shape);
#if ROOT_VERSION_CODE > ROOT_VERSION(6,21,0)
else if (isa == TGeoTessellated::Class())
solid = convertShape<TGeoTessellated>(shape);
else if (isa == TGeoScaledShape::Class()) {
TGeoScaledShape* sh = (TGeoScaledShape*) shape;
#if G4VERSION_NUMBER >= 1030
const double* vals = sh->GetScale()->GetScale();
G4Scale3D scal(vals[0], vals[1], vals[2]);
G4VSolid* g4solid = (G4VSolid*)handleSolid(sol->GetName(), sol);
if ( scal.xx()>0e0 && scal.yy()>0e0 && scal.zz()>0e0 )
solid = new G4ScaledSolid(sh->GetName(), g4solid, scal);
else
solid = new G4ReflectedSolid(g4solid->GetName()+"_refl", g4solid, scal);
#else
except("Geant4Converter","++ TGeoScaledShape are only supported by Geant4 for versions >= 10.3");
#endif
}
else if (isa == TGeoCompositeShape::Class()) {
const TGeoCompositeShape* sh = (const TGeoCompositeShape*) shape;
const TGeoBoolNode* boolean = sh->GetBoolNode();
TGeoBoolNode::EGeoBoolType oper = boolean->GetBooleanOperator();
TGeoMatrix* matrix = boolean->GetRightMatrix();
G4VSolid* left = (G4VSolid*) handleSolid(name + "_left", boolean->GetLeftShape());
G4VSolid* right = (G4VSolid*) handleSolid(name + "_right", boolean->GetRightShape());
except("Geant4Converter","++ No left Geant4 Solid present for composite shape: %s",name.c_str());
except("Geant4Converter","++ No right Geant4 Solid present for composite shape: %s",name.c_str());
//specific case!
//Ellipsoid tag preparing
//if left == TGeoScaledShape AND right == TGeoBBox
// AND if TGeoScaledShape->GetShape == TGeoSphere
TGeoShape* ls = boolean->GetLeftShape();
TGeoShape* rs = boolean->GetRightShape();
if (strcmp(ls->ClassName(), "TGeoScaledShape") == 0 &&
Markus Frank
committed
strcmp(rs->ClassName(), "TGeoBBox") == 0) {
if (strcmp(((TGeoScaledShape *)ls)->GetShape()->ClassName(), "TGeoSphere") == 0) {
if (oper == TGeoBoolNode::kGeoIntersection) {
TGeoScaledShape* lls = (TGeoScaledShape *)ls;
TGeoBBox* rrs = (TGeoBBox*)rs;
Markus Frank
committed
double sx = lls->GetScale()->GetScale()[0];
double sy = lls->GetScale()->GetScale()[1];
Markus Frank
committed
double radius = ((TGeoSphere *)lls->GetShape())->GetRmax();
Markus Frank
committed
double dz = rrs->GetDZ();
double zorig = rrs->GetOrigin()[2];
double zcut2 = dz + zorig;
double zcut1 = 2 * zorig - zcut2;
Markus Frank
committed
solid = new G4Ellipsoid(name,
sx * radius * CM_2_MM,
sy * radius * CM_2_MM,
radius * CM_2_MM,
zcut1 * CM_2_MM,
zcut2 * CM_2_MM);
data().g4Solids[shape] = solid;
return solid;
}
}
if (matrix->IsRotation()) {
MyTransform3D transform(matrix->GetTranslation(),matrix->GetRotationMatrix());
if (oper == TGeoBoolNode::kGeoSubtraction)
solid = new G4SubtractionSolid(name, left, right, transform);
else if (oper == TGeoBoolNode::kGeoUnion)
solid = new G4UnionSolid(name, left, right, transform);
else if (oper == TGeoBoolNode::kGeoIntersection)
solid = new G4IntersectionSolid(name, left, right, transform);
const Double_t *t = matrix->GetTranslation();
G4ThreeVector transform(t[0] * CM_2_MM, t[1] * CM_2_MM, t[2] * CM_2_MM);
if (oper == TGeoBoolNode::kGeoSubtraction)
solid = new G4SubtractionSolid(name, left, right, 0, transform);
else if (oper == TGeoBoolNode::kGeoUnion)
solid = new G4UnionSolid(name, left, right, 0, transform);
else if (oper == TGeoBoolNode::kGeoIntersection)
solid = new G4IntersectionSolid(name, left, right, 0, transform);
if (!solid)
except("Geant4Converter","++ Failed to handle unknown solid shape: %s of type %s",
name.c_str(), isa->GetName());
printout(lvl,"Geant4Converter","++ Successessfully converted shape [%p] of type:%s to %s.",
solid,isa->GetName(),typeName(typeid(*solid)).c_str());
data().g4Solids[shape] = solid;
}
return solid;
}
/// Dump logical volume in GDML format to output stream
void* Geant4Converter::handleVolume(const string& name, const TGeoVolume* volume) const {
Geant4GeometryInfo& info = data();
PrintLevel lvl = debugVolumes ? ALWAYS : outputLevel;
Markus Frank
committed
Geant4GeometryMaps::VolumeMap::const_iterator volIt = info.g4Volumes.find(volume);
Volume _v(volume);
if ( _v.testFlagBit(Volume::VETO_SIMU) ) {
printout(lvl, "Geant4Converter", "++ Volume %s not converted [Veto'ed for simulation]",volume->GetName());
return 0;
}
else if (volIt == info.g4Volumes.end() ) {
string n = v->GetName();
TGeoMedium* med = v->GetMedium();
TGeoShape* sh = v->GetShape();
G4VSolid* solid = (G4VSolid*) handleSolid(sh->GetName(), sh);
G4Material* medium = 0;
bool is_assembly = sh->IsA() == TGeoShapeAssembly::Class() || v->IsA() == TGeoVolumeAssembly::Class();
LimitSet lim = _v.limitSet();
Markus Frank
committed
G4UserLimits* user_limits = 0;
if (!user_limits) {
throw runtime_error("G4Cnv::volume[" + name + "]: + FATAL Failed to "
Markus Frank
committed
"access Geant4 user limits.");
VisAttr vis = _v.visAttributes();
G4VisAttributes* vis_attr = 0;
vis_attr = (G4VisAttributes*)handleVis(vis.name(), vis);
Markus Frank
committed
}
Markus Frank
committed
G4Region* region = 0;
if (!region) {
throw runtime_error("G4Cnv::volume[" + name + "]: + FATAL Failed to "
Markus Frank
committed
"access Geant4 region.");
Markus Frank
committed
}
Markus Frank
committed
printout(lvl, "Geant4Converter", "++ Convert Volume %-32s: %p %s/%s assembly:%s",
n.c_str(), v, sh->IsA()->GetName(), v->IsA()->GetName(), yes_no(is_assembly));
//info.g4AssemblyVolumes[v] = new Geant4AssemblyVolume();
Markus Frank
committed
}
medium = (G4Material*) handleMaterial(med->GetName(), Material(med));
if (!solid) {
throw runtime_error("G4Converter: No Geant4 Solid present for volume:" + n);
Markus Frank
committed
}
if (!medium) {
throw runtime_error("G4Converter: No Geant4 material present for volume:" + n);
Markus Frank
committed
}
printout(lvl, "Geant4Converter", "++ Volume + Apply LIMITS settings:%-24s to volume %s.",
Markus Frank
committed
lim.name(), _v.name());
Markus Frank
committed
}
Markus Frank
committed
G4LogicalVolume* vol = new G4LogicalVolume(solid, medium, n, 0, 0, user_limits);
printout(lvl, "Geant4Converter", "++ Volume + Apply REGION settings: %s to volume %s.",
Markus Frank
committed
reg.name(), _v.name());
Markus Frank
committed
vol->SetRegion(region);
region->AddRootLogicalVolume(vol);
Markus Frank
committed
vol->SetVisAttributes(vis_attr);
}
info.g4Volumes[v] = vol;
printout(lvl, "Geant4Converter", "++ Volume + %s converted: %p ---> G4: %p", n.c_str(), v, vol);
return nullptr;
/// Dump logical volume in GDML format to output stream
void* Geant4Converter::collectVolume(const string& /* name */, const TGeoVolume* volume) const {
Geant4GeometryInfo& info = data();
Region reg = _v.region();
LimitSet lim = _v.limitSet();
SensitiveDetector det = _v.sensitiveDetector();
/// Dump volume placement in GDML format to output stream
void* Geant4Converter::handleAssembly(const string& name, const TGeoNode* node) const {
TGeoVolume* mot_vol = node->GetVolume();
PrintLevel lvl = debugVolumes ? ALWAYS : outputLevel;
if ( mot_vol->IsA() != TGeoVolumeAssembly::Class() ) {
Volume _v(mot_vol);
if ( _v.testFlagBit(Volume::VETO_SIMU) ) {
printout(lvl, "Geant4Converter", "++ AssemblyNode %s not converted [Veto'ed for simulation]",node->GetName());
return 0;
}
Geant4GeometryInfo& info = data();
Geant4AssemblyVolume* g4 = info.g4AssemblyVolumes[node];
if ( g4 ) {
printout(ALWAYS, "Geant4Converter", "+++ Assembly: **** : Re-using existing assembly: %s",node->GetName());
}
if ( !g4 ) {
g4 = new Geant4AssemblyVolume();
for(Int_t i=0; i < mot_vol->GetNdaughters(); ++i) {
TGeoNode* dau = mot_vol->GetNode(i);
TGeoVolume* dau_vol = dau->GetVolume();
TGeoMatrix* tr = dau->GetMatrix();
MyTransform3D transform(tr->GetTranslation(),tr->IsRotation() ? tr->GetRotationMatrix() : s_identity_rot);
if ( is_left_handed(tr) ) {
G4Scale3D scale;
G4Rotate3D rot;
G4Translate3D trans;
transform.getDecomposition(scale, rot, trans);
printout(debugReflections ? ALWAYS : lvl, "Geant4Converter",
"+++ Placing reflected ASSEMBLY. dau:%s to mother %s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
dau_vol->GetName(), mot_vol->GetName(),
transform.dx(), transform.dy(), transform.dz(),
scale.xx(), scale.yy(), scale.zz());
if ( dau_vol->IsA() == TGeoVolumeAssembly::Class() ) {
Geant4GeometryMaps::AssemblyMap::iterator ia = info.g4AssemblyVolumes.find(dau);
if ( ia == info.g4AssemblyVolumes.end() ) {
Markus Frank
committed
printout(FATAL, "Geant4Converter", "+++ Invalid child assembly at %s : %d parent: %s child:%s",
__FILE__, __LINE__, name.c_str(), dau->GetName());
Markus Frank
committed
}
g4->placeAssembly(dau,(*ia).second,transform);
printout(lvl, "Geant4Converter", "+++ Assembly: AddPlacedAssembly : dau:%s "
Markus Frank
committed
"to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f",
dau_vol->GetName(), mot_vol->GetName(),
transform.dx(), transform.dy(), transform.dz());
Geant4GeometryMaps::VolumeMap::iterator iv = info.g4Volumes.find(dau_vol);
if ( iv == info.g4Volumes.end() ) {
printout(FATAL,"Geant4Converter", "+++ Invalid child volume at %s : %d parent: %s child:%s",
__FILE__, __LINE__, name.c_str(), dau->GetName());
except("Geant4Converter", "+++ Invalid child volume at %s : %d parent: %s child:%s",
__FILE__, __LINE__, name.c_str(), dau->GetName());
Markus Frank
committed
}
g4->placeVolume(dau,(*iv).second, transform);
printout(lvl, "Geant4Converter", "+++ Assembly: AddPlacedVolume : dau:%s "
Markus Frank
committed
"to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f",
dau_vol->GetName(), mot_vol->GetName(),
transform.dx(), transform.dy(), transform.dz());
info.g4AssemblyVolumes[node] = g4;
/// Dump volume placement in GDML format to output stream
void* Geant4Converter::handlePlacement(const string& name, const TGeoNode* node) const {
Geant4GeometryInfo& info = data();
PrintLevel lvl = debugPlacements ? ALWAYS : outputLevel;
Markus Frank
committed
Geant4GeometryMaps::PlacementMap::const_iterator g4it = info.g4Placements.find(node);
G4VPhysicalVolume* g4 = (g4it == info.g4Placements.end()) ? 0 : (*g4it).second;
TGeoVolume* vol = node->GetVolume();
Volume _v(vol);
if ( _v.testFlagBit(Volume::VETO_SIMU) ) {
printout(lvl, "Geant4Converter", "++ Placement %s not converted [Veto'ed for simulation]",node->GetName());
return 0;
}
if (!g4) {
TGeoVolume* mot_vol = node->GetMotherVolume();
TGeoMatrix* tr = node->GetMatrix();
if (!tr) {
except("Geant4Converter",
"++ Attempt to handle placement without transformation:%p %s of type %s vol:%p",
node, node->GetName(), node->IsA()->GetName(), vol);
except("Geant4Converter", "++ Unknown G4 volume:%p %s of type %s ptr:%p",
node, node->GetName(), node->IsA()->GetName(), vol);
int copy = node->GetNumber();
bool node_is_reflected = is_left_handed(tr);
bool node_is_assembly = vol->IsA() == TGeoVolumeAssembly::Class();
bool mother_is_assembly = mot_vol ? mot_vol->IsA() == TGeoVolumeAssembly::Class() : false;
MyTransform3D transform(tr->GetTranslation(),tr->IsRotation() ? tr->GetRotationMatrix() : s_identity_rot);
Markus Frank
committed
Geant4GeometryMaps::VolumeMap::const_iterator volIt = info.g4Volumes.find(mot_vol);
if ( mother_is_assembly ) {
Markus Frank
committed
//
// Mother is an assembly:
// Nothing to do here, because:
// -- placed volumes were already added before in "handleAssembly"
// -- imprint cannot be made, because this requires a logical volume as a mother
//
printout(lvl, "Geant4Converter", "+++ Assembly: **** : dau:%s "
Markus Frank
committed
"to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f",
vol->GetName(), mot_vol->GetName(),
transform.dx(), transform.dy(), transform.dz());
Markus Frank
committed
}
G4Scale3D scale;
G4Rotate3D rot;
G4Translate3D trans;
if ( node_is_assembly ) {
Markus Frank
committed
//
// Node is an assembly:
// Imprint the assembly. The mother MUST already be transformed.
//
transform.getDecomposition(scale, rot, trans);
printout(lvl, "Geant4Converter", "+++ Assembly: makeImprint: dau:%-12s %s in mother %-12s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
node->GetName(), node_is_reflected ? "(REFLECTED)" : "",
mot_vol ? mot_vol->GetName() : "<unknown>",
transform.dx(), transform.dy(), transform.dz(),
scale.xx(), scale.yy(), scale.zz());
Markus Frank
committed
Geant4AssemblyVolume* ass = (Geant4AssemblyVolume*)info.g4AssemblyVolumes[node];
Geant4AssemblyVolume::Chain chain;
chain.emplace_back(node);
Markus Frank
committed
ass->imprint(info,node,chain,ass,(*volIt).second, transform, copy, checkOverlaps);
Markus Frank
committed
return 0;
else if ( node != info.manager->GetTopNode() && volIt == info.g4Volumes.end() ) {
Markus Frank
committed
throw logic_error("Geant4Converter: Invalid mother volume found!");
G4LogicalVolume* g4vol = info.g4Volumes[vol];
G4LogicalVolume* g4mot = info.g4Volumes[mot_vol];
G4PhysicalVolumesPair pvPlaced =
G4ReflectionFactory::Instance()->Place(transform, // no rotation
name, // its name
g4vol, // its logical volume
g4mot, // its mother (logical) volume
false, // no boolean operations
copy, // its copy number
checkOverlaps);
transform.getDecomposition(scale, rot, trans);
printout(debugReflections ? ALWAYS : lvl, "Geant4Converter",
"+++ Place %svolume %-12s in mother %-12s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
node_is_reflected ? "REFLECTED " : "", _v.name(),
mot_vol ? mot_vol->GetName() : "<unknown>",
transform.dx(), transform.dy(), transform.dz(),
scale.xx(), scale.yy(), scale.zz());
// First 2 cases can be combined.
// Leave them separated for debugging G4ReflectionFactory for now...
if ( node_is_reflected && !pvPlaced.second )
return info.g4Placements[node] = pvPlaced.first;
else if ( !node_is_reflected && !pvPlaced.second )
return info.g4Placements[node] = pvPlaced.first;
//G4LogicalVolume* g4refMoth = G4ReflectionFactory::Instance()->GetReflectedLV(g4mot);
// Now deal with valid pvPlaced.second ...
if ( node_is_reflected )
return info.g4Placements[node] = pvPlaced.first;
else if ( !node_is_reflected )
return info.g4Placements[node] = pvPlaced.first;
g4 = pvPlaced.second ? pvPlaced.second : pvPlaced.first;
info.g4Placements[node] = g4;
printout(ERROR, "Geant4Converter", "++ DEAD code. Should not end up here!");
printout(ERROR, "Geant4Converter", "++ Attempt to DOUBLE-place physical volume: %s No:%d", name.c_str(), node->GetNumber());
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
void* Geant4Converter::handlePlacement2(const std::pair<const TGeoNode*,const TGeoNode*>& n) const {
const TGeoNode* par_node = n.first;
const TGeoNode* node = n.second;
Geant4GeometryInfo& info = data();
PrintLevel lvl = debugPlacements ? ALWAYS : outputLevel;
Geant4GeometryMaps::PlacementMap::const_iterator g4it = info.g4Placements.find(node);
G4VPhysicalVolume* g4 = (g4it == info.g4Placements.end()) ? 0 : (*g4it).second;
TGeoVolume* vol = node->GetVolume();
Volume _v(vol);
if ( _v.testFlagBit(Volume::VETO_SIMU) ) {
printout(lvl, "Geant4Converter", "++ Placement %s not converted [Veto'ed for simulation]",node->GetName());
return 0;
}
g4 = nullptr;
if (!g4) {
TGeoVolume* mot_vol = par_node ? par_node->GetVolume() : node->GetMotherVolume();
TGeoMatrix* tr = node->GetMatrix();
if (!tr) {
except("Geant4Converter",
"++ Attempt to handle placement without transformation:%p %s of type %s vol:%p",
node, node->GetName(), node->IsA()->GetName(), vol);
}
else if (0 == vol) {
except("Geant4Converter", "++ Unknown G4 volume:%p %s of type %s ptr:%p",
node, node->GetName(), node->IsA()->GetName(), vol);
}
else {
int copy = node->GetNumber();
bool node_is_reflected = is_left_handed(tr);
bool node_is_assembly = vol->IsA() == TGeoVolumeAssembly::Class();
bool mother_is_assembly = mot_vol ? mot_vol->IsA() == TGeoVolumeAssembly::Class() : false;
MyTransform3D transform(tr->GetTranslation(),tr->IsRotation() ? tr->GetRotationMatrix() : s_identity_rot);
Geant4GeometryMaps::VolumeMap::const_iterator volIt = info.g4Volumes.find(mot_vol);
if ( mother_is_assembly ) {
//
// Mother is an assembly:
// Nothing to do here, because:
// -- placed volumes were already added before in "handleAssembly"
// -- imprint cannot be made, because this requires a logical volume as a mother
//
printout(lvl, "Geant4Converter", "+++ Assembly: **** : dau:%s "
"to mother %s Tr:x=%8.3f y=%8.3f z=%8.3f",
vol->GetName(), mot_vol->GetName(),
transform.dx(), transform.dy(), transform.dz());
return 0;
}
G4Scale3D scale;
G4Rotate3D rot;
G4Translate3D trans;
if ( node_is_assembly ) {
//
// Node is an assembly:
// Imprint the assembly. The mother MUST already be transformed.
//
transform.getDecomposition(scale, rot, trans);
printout(lvl, "Geant4Converter", "+++ Assembly: makeImprint: dau:%-12s %s in mother %-12s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
node->GetName(), node_is_reflected ? "(REFLECTED)" : "",
mot_vol ? mot_vol->GetName() : "<unknown>",
transform.dx(), transform.dy(), transform.dz(),
scale.xx(), scale.yy(), scale.zz());
Geant4AssemblyVolume* ass = (Geant4AssemblyVolume*)info.g4AssemblyVolumes[node];
Geant4AssemblyVolume::Chain chain;
chain.emplace_back(node);
ass->imprint(info,node,chain,ass,(*volIt).second, transform, copy, checkOverlaps);
return 0;
}
else if ( node != info.manager->GetTopNode() && volIt == info.g4Volumes.end() ) {
throw logic_error("Geant4Converter: Invalid mother volume found!");
}
G4LogicalVolume* g4vol = info.g4Volumes[vol];
G4LogicalVolume* g4mot = info.g4Volumes[mot_vol];
G4PhysicalVolumesPair pvPlaced =
G4ReflectionFactory::Instance()->Place(transform, // no rotation
node->GetName(), // its name
g4vol, // its logical volume
g4mot, // its mother (logical) volume
false, // no boolean operations
copy, // its copy number
checkOverlaps);
transform.getDecomposition(scale, rot, trans);
printout(debugReflections ? ALWAYS : lvl, "Geant4Converter",
"+++ Place %svolume %-12s in mother %-12s "
"Tr:x=%8.1f y=%8.1f z=%8.1f Scale:x=%4.2f y=%4.2f z=%4.2f",
node_is_reflected ? "REFLECTED " : "", _v.name(),
mot_vol ? mot_vol->GetName() : "<unknown>",
transform.dx(), transform.dy(), transform.dz(),
scale.xx(), scale.yy(), scale.zz());
// First 2 cases can be combined.
// Leave them separated for debugging G4ReflectionFactory for now...
if ( node_is_reflected && !pvPlaced.second )
return info.g4Placements[node] = pvPlaced.first;
else if ( !node_is_reflected && !pvPlaced.second )
return info.g4Placements[node] = pvPlaced.first;
//G4LogicalVolume* g4refMoth = G4ReflectionFactory::Instance()->GetReflectedLV(g4mot);
// Now deal with valid pvPlaced.second ...
if ( node_is_reflected )
return info.g4Placements[node] = pvPlaced.first;
else if ( !node_is_reflected )
return info.g4Placements[node] = pvPlaced.first;
g4 = pvPlaced.second ? pvPlaced.second : pvPlaced.first;
}
info.g4Placements[node] = g4;
printout(ERROR, "Geant4Converter", "++ DEAD code. Should not end up here!");
}
else {
printout(ERROR, "Geant4Converter", "++ Attempt to DOUBLE-place physical volume: %s No:%d",
node->GetName(), node->GetNumber());
}
return g4;
}
/// Convert the geometry type region into the corresponding Geant4 object(s).
void* Geant4Converter::handleRegion(Region region, const set<const TGeoVolume*>& /* volumes */) const {
G4Region* g4 = data().g4Regions[region];
Markus Frank
committed
PrintLevel lvl = debugRegions ? ALWAYS : outputLevel;
Region r = region;
Markus Frank
committed
g4 = new G4Region(r.name());
// create region info with storeSecondaries flag
Andre Sailer
committed
if( not r.wasThresholdSet() and r.storeSecondaries() ) {
throw runtime_error("G4Region: StoreSecondaries is True, but no explicit threshold set:");
}
printout(lvl, "Geant4Converter", "++ Setting up region: %s", r.name());
G4UserRegionInformation* info = new G4UserRegionInformation();
info->region = r;
info->threshold = r.threshold()*CLHEP::MeV/units::MeV;
info->storeSecondaries = r.storeSecondaries();
g4->SetUserInformation(info);
printout(lvl, "Geant4Converter", "++ Converted region settings of:%s.", r.name());
vector < string > &limits = r.limits();
G4ProductionCuts* cuts = 0;
// set production cut
if( not r.useDefaultCut() ) {
cuts = new G4ProductionCuts();
cuts->SetProductionCut(r.cut()*CLHEP::mm/units::mm);
printout(lvl, "Geant4Converter", "++ %s: Using default cut: %f [mm]",
r.name(), r.cut()*CLHEP::mm/units::mm);
}
for (const auto& nam : limits ) {
LimitSet ls = m_detDesc.limitSet(nam);
const LimitSet::Set& cts = ls.cuts();
for (const auto& c : cts ) {
if ( c.particles == "*" ) pid = -1;
else if ( c.particles == "e-" ) pid = idxG4ElectronCut;
else if ( c.particles == "e+" ) pid = idxG4PositronCut;
else if ( c.particles == "e[+-]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
else if ( c.particles == "e[-+]" ) pid = -idxG4PositronCut-idxG4ElectronCut;
else if ( c.particles == "gamma" ) pid = idxG4GammaCut;
else if ( c.particles == "proton" ) pid = idxG4ProtonCut;
else throw runtime_error("G4Region: Invalid production cut particle-type:" + c.particles);
if ( !cuts ) cuts = new G4ProductionCuts();
if ( pid == -(idxG4PositronCut+idxG4ElectronCut) ) {
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4PositronCut);
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, idxG4ElectronCut);
}
else {
cuts->SetProductionCut(c.value*CLHEP::mm/units::mm, pid);
}
printout(lvl, "Geant4Converter", "++ %s: Set cut [%s/%d] = %f [mm]",
r.name(), c.particles.c_str(), pid, c.value*CLHEP::mm/units::mm);
}
const auto& lm = data().g4Limits;
for (const auto& j : lm ) {
if (nam == j.first->GetName()) {
g4->SetUserLimits(j.second);
printout(lvl, "Geant4Converter", "++ %s: Set limits %s to region type %s",
r.name(), nam.c_str(), j.second->GetType().c_str());
found = true;
break;
}
}
if ( found ) {
throw runtime_error("G4Region: Failed to resolve user limitset:" + nam);
/// Assign cuts to region if they were created
if ( cuts ) g4->SetProductionCuts(cuts);
data().g4Regions[region] = g4;
}
return g4;
}
/// Convert the geometry type LimitSet into the corresponding Geant4 object(s).
void* Geant4Converter::handleLimitSet(LimitSet limitset, const set<const TGeoVolume*>& /* volumes */) const {
G4UserLimits* g4 = data().g4Limits[limitset];
if (!g4) {
1155
1156
1157
1158
1159
1160
1161
1162
1163
1164
1165
1166
1167
1168
1169
1170
1171
1172
1173
1174
1175
1176
1177
1178
1179
1180
1181
1182
struct LimitPrint {
const LimitSet& ls;
LimitPrint(const LimitSet& lset) : ls(lset) {}
const LimitPrint& operator()(const std::string& pref, const Geant4UserLimits::Handler& h) const {
if ( !h.particleLimits.empty() ) {
printout(ALWAYS,"Geant4Converter",
"+++ LimitSet: Limit %s.%s applied for particles:",ls.name(), pref.c_str());
for(const auto& p : h.particleLimits)
printout(ALWAYS,"Geant4Converter","+++ LimitSet: Particle type: %-18s PDG: %-6d : %f",
p.first->GetParticleName().c_str(), p.first->GetPDGEncoding(), p.second);
}
else {
printout(ALWAYS,"Geant4Converter",
"+++ LimitSet: Limit %s.%s NOT APPLIED.",ls.name(), pref.c_str());
}
return *this;
}
};
Geant4UserLimits* limits = new Geant4UserLimits(limitset);
g4 = limits;
if ( debugRegions ) {
LimitPrint print(limitset);
print("maxTime", limits->maxTime)
("minEKine", limits->minEKine)
("minRange", limits->minRange)
("maxStepLength", limits->maxStepLength)
("maxTrackLength",limits->maxTrackLength);
}
data().g4Limits[limitset] = g4;
}
return g4;
}
/// Convert the geometry visualisation attributes to the corresponding Geant4 object(s).
void* Geant4Converter::handleVis(const string& /* name */, VisAttr attr) const {
Geant4GeometryInfo& info = data();
G4VisAttributes* g4 = info.g4Vis[attr];
if ( !g4 ) {
int style = attr.lineStyle();
attr.rgb(red, green, blue);
g4 = new G4VisAttributes(attr.visible(), G4Colour(red, green, blue, attr.alpha()));
//g4->SetLineWidth(attr->GetLineWidth());
g4->SetDaughtersInvisible(!attr.showDaughters());
if ( style == VisAttr::SOLID ) {
g4->SetLineStyle(G4VisAttributes::unbroken);
g4->SetForceWireframe(false);
g4->SetForceSolid(true);
}
else if ( style == VisAttr::WIREFRAME || style == VisAttr::DASHED ) {
g4->SetLineStyle(G4VisAttributes::dashed);
g4->SetForceSolid(false);
g4->SetForceWireframe(true);
}
}
return g4;
}
/// Handle the geant 4 specific properties
void Geant4Converter::handleProperties(Detector::Properties& prp) const {
map < string, string > processors;
static int s_idd = 9999999;
string id;
for (Detector::Properties::const_iterator i = prp.begin(); i != prp.end(); ++i) {
const string& nam = (*i).first;
const Detector::PropertyValues& vals = (*i).second;
if (nam.substr(0, 6) == "geant4") {
Detector::PropertyValues::const_iterator id_it = vals.find("id");
if (id_it != vals.end()) {
id = (*id_it).second;
::snprintf(txt, sizeof(txt), "%d", ++s_idd);
processors.emplace(id, nam);
for (map<string, string>::const_iterator i = processors.begin(); i != processors.end(); ++i) {
const Detector::PropertyValues& vals = prp[nam];
string type = vals.find("type")->second;
string tag = type + "_Geant4_action";
Detector* detPtr = const_cast<Detector*>(&m_detDesc);
long result = PluginService::Create<long>(tag, detPtr, hdlr, &vals);
throw runtime_error("Failed to locate plugin to interprete files of type"
Markus Frank
committed
" \"" + tag + "\" - no factory:" + type);
result = *(long*) result;
if (result != 1) {
throw runtime_error("Failed to invoke the plugin " + tag + " of type " + type);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "+++++ Executed Successfully Geant4 setup module *%s*.", type.c_str());
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
/// Convert the geometry type material into the corresponding Geant4 object(s).
void* Geant4Converter::handleMaterialProperties(TObject* matrix) const {
Geant4GeometryInfo& info = data();
Geant4GeometryInfo::PropertyVector* g4 = info.g4OpticalProperties[gdmlMat];
Markus Frank
committed
PrintLevel lvl = debugMaterials ? ALWAYS : outputLevel;
g4 = new Geant4GeometryInfo::PropertyVector();
size_t rows = gdmlMat->GetRows();
g4->name = gdmlMat->GetName();
g4->title = gdmlMat->GetTitle();
g4->bins.reserve(rows);
g4->values.reserve(rows);
for(size_t i=0; i<rows; ++i) {
g4->bins.emplace_back(gdmlMat->Get(i,0) /* *CLHEP::eV/units::eV */);
g4->values.emplace_back(gdmlMat->Get(i,1));
Markus Frank
committed
printout(lvl, "Geant4Converter", "++ Successfully converted material property:%s : %s [%ld rows]",
gdmlMat->GetName(), gdmlMat->GetTitle(), rows);
info.g4OpticalProperties[gdmlMat] = g4;
1276
1277
1278
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292
1293
1294
1295
1296
1297
1298
1299
1300
1301
1302
1303
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325
1326
1327
1328
1329
1330
1331
1332
1333
1334
1335
1336
1337
1338
1339
1340
1341
1342
1343
1344
1345
1346
}
return g4;
}
static G4OpticalSurfaceFinish geant4_surface_finish(TGeoOpticalSurface::ESurfaceFinish f) {
#define TO_G4_FINISH(x) case TGeoOpticalSurface::kF##x : return x;
switch(f) {
TO_G4_FINISH(polished); // smooth perfectly polished surface
TO_G4_FINISH(polishedfrontpainted); // smooth top-layer (front) paint
TO_G4_FINISH(polishedbackpainted); // same is 'polished' but with a back-paint
TO_G4_FINISH(ground); // rough surface
TO_G4_FINISH(groundfrontpainted); // rough top-layer (front) paint
TO_G4_FINISH(groundbackpainted); // same as 'ground' but with a back-paint
TO_G4_FINISH(polishedlumirrorair); // mechanically polished surface, with lumirror
TO_G4_FINISH(polishedlumirrorglue); // mechanically polished surface, with lumirror & meltmount
TO_G4_FINISH(polishedair); // mechanically polished surface
TO_G4_FINISH(polishedteflonair); // mechanically polished surface, with teflon
TO_G4_FINISH(polishedtioair); // mechanically polished surface, with tio paint
TO_G4_FINISH(polishedtyvekair); // mechanically polished surface, with tyvek
TO_G4_FINISH(polishedvm2000air); // mechanically polished surface, with esr film
TO_G4_FINISH(polishedvm2000glue); // mechanically polished surface, with esr film & meltmount
TO_G4_FINISH(etchedlumirrorair); // chemically etched surface, with lumirror
TO_G4_FINISH(etchedlumirrorglue); // chemically etched surface, with lumirror & meltmount
TO_G4_FINISH(etchedair); // chemically etched surface
TO_G4_FINISH(etchedteflonair); // chemically etched surface, with teflon
TO_G4_FINISH(etchedtioair); // chemically etched surface, with tio paint
TO_G4_FINISH(etchedtyvekair); // chemically etched surface, with tyvek
TO_G4_FINISH(etchedvm2000air); // chemically etched surface, with esr film
TO_G4_FINISH(etchedvm2000glue); // chemically etched surface, with esr film & meltmount
TO_G4_FINISH(groundlumirrorair); // rough-cut surface, with lumirror
TO_G4_FINISH(groundlumirrorglue); // rough-cut surface, with lumirror & meltmount
TO_G4_FINISH(groundair); // rough-cut surface
TO_G4_FINISH(groundteflonair); // rough-cut surface, with teflon
TO_G4_FINISH(groundtioair); // rough-cut surface, with tio paint
TO_G4_FINISH(groundtyvekair); // rough-cut surface, with tyvek
TO_G4_FINISH(groundvm2000air); // rough-cut surface, with esr film
TO_G4_FINISH(groundvm2000glue); // rough-cut surface, with esr film & meltmount
// for DAVIS model
TO_G4_FINISH(Rough_LUT); // rough surface
TO_G4_FINISH(RoughTeflon_LUT); // rough surface wrapped in Teflon tape
TO_G4_FINISH(RoughESR_LUT); // rough surface wrapped with ESR
TO_G4_FINISH(RoughESRGrease_LUT); // rough surface wrapped with ESR and coupled with opical grease
TO_G4_FINISH(Polished_LUT); // polished surface
TO_G4_FINISH(PolishedTeflon_LUT); // polished surface wrapped in Teflon tape
TO_G4_FINISH(PolishedESR_LUT); // polished surface wrapped with ESR
TO_G4_FINISH(PolishedESRGrease_LUT); // polished surface wrapped with ESR and coupled with opical grease
TO_G4_FINISH(Detector_LUT); // polished surface with optical grease
default:
printout(ERROR,"Geant4Surfaces","++ Unknown finish style: %d [%s]. Assume polished!",
int(f), TGeoOpticalSurface::FinishToString(f));
return polished;
}
#undef TO_G4_FINISH
}
static G4SurfaceType geant4_surface_type(TGeoOpticalSurface::ESurfaceType t) {
#define TO_G4_TYPE(x) case TGeoOpticalSurface::kT##x : return x;
switch(t) {
TO_G4_TYPE(dielectric_metal); // dielectric-metal interface
TO_G4_TYPE(dielectric_dielectric); // dielectric-dielectric interface
TO_G4_TYPE(dielectric_LUT); // dielectric-Look-Up-Table interface
TO_G4_TYPE(dielectric_LUTDAVIS); // dielectric-Look-Up-Table DAVIS interface
TO_G4_TYPE(dielectric_dichroic); // dichroic filter interface
TO_G4_TYPE(firsov); // for Firsov Process
TO_G4_TYPE(x_ray); // for x-ray mirror process
default:
printout(ERROR,"Geant4Surfaces","++ Unknown surface type: %d [%s]. Assume dielectric_metal!",
int(t), TGeoOpticalSurface::TypeToString(t));
return dielectric_metal;
}
#undef TO_G4_TYPE
}
static G4OpticalSurfaceModel geant4_surface_model(TGeoOpticalSurface::ESurfaceModel surfMod) {
#define TO_G4_MODEL(x) case TGeoOpticalSurface::kM##x : return x;
TO_G4_MODEL(glisur); // original GEANT3 model
TO_G4_MODEL(unified); // UNIFIED model
TO_G4_MODEL(LUT); // Look-Up-Table model
TO_G4_MODEL(DAVIS); // DAVIS model
TO_G4_MODEL(dichroic); // dichroic filter
default:
printout(ERROR,"Geant4Surfaces","++ Unknown surface model: %d [%s]. Assume glisur!",
int(surfMod), TGeoOpticalSurface::ModelToString(surfMod));
return glisur;
}
#undef TO_G4_MODEL
}
/// Convert the optical surface to Geant4
void* Geant4Converter::handleOpticalSurface(TObject* surface) const {
TGeoOpticalSurface* optSurf = (TGeoOpticalSurface*)surface;
Geant4GeometryInfo& info = data();
G4OpticalSurface* g4 = info.g4OpticalSurfaces[optSurf];
G4SurfaceType type = geant4_surface_type(optSurf->GetType());
G4OpticalSurfaceModel model = geant4_surface_model(optSurf->GetModel());
G4OpticalSurfaceFinish finish = geant4_surface_finish(optSurf->GetFinish());
g4 = new G4OpticalSurface(optSurf->GetName(), model, finish, type, optSurf->GetValue());
g4->SetSigmaAlpha(optSurf->GetSigmaAlpha());
// not implemented: g4->SetPolish(s->GetPolish());
Markus Frank
committed
printout(debugSurfaces ? ALWAYS : DEBUG, "Geant4Converter",
"++ Created OpticalSurface: %-18s type:%s model:%s finish:%s",
optSurf->GetName(),
TGeoOpticalSurface::TypeToString(optSurf->GetType()),
TGeoOpticalSurface::ModelToString(optSurf->GetModel()),
TGeoOpticalSurface::FinishToString(optSurf->GetFinish()));
G4MaterialPropertiesTable* tab = 0;
Markus Frank
committed
for(TObject* obj = it.Next(); obj; obj = it.Next()) {
TNamed* named = (TNamed*)obj;
TGDMLMatrix* matrix = info.manager->GetGDMLMatrix(named->GetTitle());
if ( 0 == tab ) {
tab = new G4MaterialPropertiesTable();
g4->SetMaterialPropertiesTable(tab);
}
Geant4GeometryInfo::PropertyVector* v =
(Geant4GeometryInfo::PropertyVector*)handleMaterialProperties(matrix);
if ( !v ) { // Error!
except("Geant4OpticalSurface","++ Failed to convert opt.surface %s. Property table %s is not defined!",
optSurf->GetName(), named->GetTitle());
}
int idx = tab->GetPropertyIndex(named->GetName(), false);
if ( idx < 0 ) {
printout(ERROR, "Geant4Converter", "++ UNKNOWN Geant4 Property: %-20s [IGNORED]", named->GetName());
continue;
// We need to convert the property from TGeo units to Geant4 units
auto conv = g4PropertyConversion(idx);
vector<double> bins(v->bins), vals(v->values);
for(size_t i=0, count=v->bins.size(); i<count; ++i)
bins[i] *= conv.first, vals[i] *= conv.second;
G4MaterialPropertyVector* vec = new G4MaterialPropertyVector(&bins[0], &vals[0], bins.size());
tab->AddProperty(named->GetName(), vec);
Markus Frank
committed
printout(debugSurfaces ? ALWAYS : DEBUG, "Geant4Converter",
"++ Property: %-20s [%ld x %ld] --> %s",
named->GetName(), matrix->GetRows(), matrix->GetCols(), named->GetTitle());
for(size_t i=0, count=v->bins.size(); i<count; ++i)
printout(debugSurfaces ? ALWAYS : DEBUG, named->GetName(),
" Geant4: %8.3g [MeV] TGeo: %8.3g [GeV] Conversion: %8.3g",
bins[i], v->bins[i], conv.first);
}
return g4;
}
/// Convert the skin surface to Geant4
void* Geant4Converter::handleSkinSurface(TObject* surface) const {
TGeoSkinSurface* surf = (TGeoSkinSurface*)surface;
Geant4GeometryInfo& info = data();
G4LogicalSkinSurface* g4 = info.g4SkinSurfaces[surf];
if (!g4) {
G4OpticalSurface* optSurf = info.g4OpticalSurfaces[OpticalSurface(surf->GetSurface())];
G4LogicalVolume* v = info.g4Volumes[surf->GetVolume()];
g4 = new G4LogicalSkinSurface(surf->GetName(), v, optSurf);
Markus Frank
committed
printout(debugSurfaces ? ALWAYS : DEBUG, "Geant4Converter",
"++ Created SkinSurface: %-18s optical:%s",
surf->GetName(), surf->GetSurface()->GetName());
info.g4SkinSurfaces[surf] = g4;
}
return g4;
}
/// Convert the border surface to Geant4
void* Geant4Converter::handleBorderSurface(TObject* surface) const {
TGeoBorderSurface* surf = (TGeoBorderSurface*)surface;
Geant4GeometryInfo& info = data();
G4LogicalBorderSurface* g4 = info.g4BorderSurfaces[surf];
if (!g4) {
G4OpticalSurface* optSurf = info.g4OpticalSurfaces[OpticalSurface(surf->GetSurface())];
G4VPhysicalVolume* n1 = info.g4Placements[surf->GetNode1()];
G4VPhysicalVolume* n2 = info.g4Placements[surf->GetNode2()];
g4 = new G4LogicalBorderSurface(surf->GetName(), n1, n2, optSurf);
Markus Frank
committed
printout(debugSurfaces ? ALWAYS : DEBUG, "Geant4Converter",
"++ Created BorderSurface: %-18s optical:%s",
surf->GetName(), surf->GetSurface()->GetName());
info.g4BorderSurfaces[surf] = g4;
}
return g4;
}
#endif
/// Convert the geometry type SensitiveDetector into the corresponding Geant4 object(s).
Markus Frank
committed
void Geant4Converter::printSensitive(SensitiveDetector sens_det, const set<const TGeoVolume*>& /* volumes */) const {
Geant4GeometryInfo& info = data();
set<const TGeoVolume*>& volset = info.sensitives[sens_det];
SensitiveDetector sd = sens_det;
Markus Frank
committed
stringstream str;
printout(INFO, "Geant4Converter", "++ SensitiveDetector: %-18s %-20s Hits:%-16s", sd.name(), ("[" + sd.type() + "]").c_str(),
Markus Frank
committed
sd.hitsCollection().c_str());
str << " | " << "Cutoff:" << setw(6) << left << sd.energyCutoff() << setw(5) << right << volset.size()
<< " volumes ";
if (sd.region().isValid())
str << " Region:" << setw(12) << left << sd.region().name();
if (sd.limits().isValid())
str << " Limits:" << setw(12) << left << sd.limits().name();
Markus Frank
committed
str << ".";
printout(INFO, "Geant4Converter", str.str().c_str());
for (const auto i : volset ) {
map<Volume, G4LogicalVolume*>::iterator v = info.g4Volumes.find(i);
Markus Frank
committed
str.str("");
str << " | " << "Volume:" << setw(24) << left << vol->GetName() << " "
<< vol->GetNoDaughters() << " daughters.";
Markus Frank
committed
printout(INFO, "Geant4Converter", str.str().c_str());
Markus Frank
committed
string printSolid(G4VSolid* sol) {
stringstream str;
if (typeid(*sol) == typeid(G4Box)) {
const G4Box* b = (G4Box*) sol;
Markus Frank
committed
str << "++ Box: x=" << b->GetXHalfLength() << " y=" << b->GetYHalfLength() << " z=" << b->GetZHalfLength();
else if (typeid(*sol) == typeid(G4Tubs)) {
const G4Tubs* t = (const G4Tubs*) sol;
str << " Tubs: Ri=" << t->GetInnerRadius() << " Ra=" << t->GetOuterRadius() << " z/2=" << t->GetZHalfLength() << " Phi="
<< t->GetStartPhiAngle() << "..." << t->GetDeltaPhiAngle();
Markus Frank
committed
return str.str();
}
/// Print G4 placement
void* Geant4Converter::printPlacement(const string& name, const TGeoNode* node) const {
Geant4GeometryInfo& info = data();
G4VPhysicalVolume* g4 = info.g4Placements[node];
G4LogicalVolume* vol = info.g4Volumes[node->GetVolume()];
G4LogicalVolume* mot = info.g4Volumes[node->GetMotherVolume()];
G4VSolid* sol = vol->GetSolid();
G4ThreeVector tr = g4->GetObjectTranslation();
G4VSensitiveDetector* sd = vol->GetSensitiveDetector();
Markus Frank
committed
stringstream str;
str << "G4Cnv::placement: + " << name << " No:" << node->GetNumber() << " Vol:" << vol->GetName() << " Solid:"
<< sol->GetName();
Markus Frank
committed
printout(outputLevel, "G4Placement", str.str().c_str());
Markus Frank
committed
str.str("");
str << " |" << " Loc: x=" << tr.x() << " y=" << tr.y() << " z=" << tr.z();
Markus Frank
committed
printout(outputLevel, "G4Placement", str.str().c_str());
printout(outputLevel, "G4Placement", printSolid(sol).c_str());
Markus Frank
committed
str.str("");
str << " |" << " Ndau:" << vol->GetNoDaughters() << " physvols." << " Mat:" << vol->GetMaterial()->GetName()
<< " Mother:" << (char*) (mot ? mot->GetName().c_str() : "---");
Markus Frank
committed
printout(outputLevel, "G4Placement", str.str().c_str());
Markus Frank
committed
str.str("");
Markus Frank
committed
printout(outputLevel, "G4Placement", str.str().c_str());
namespace {
template <typename O, typename C, typename F> void handleRefs(const O* o, const C& c, F pmf) {
for (typename C::const_iterator i = c.begin(); i != c.end(); ++i) {
//(o->*pmf)((*i)->GetName(), *i);
(o->*pmf)("", *i);
}
template <typename O, typename C, typename F> void handle(const O* o, const C& c, F pmf) {
for (typename C::const_iterator i = c.begin(); i != c.end(); ++i) {
(o->*pmf)((*i)->GetName(), *i);
}
template <typename O, typename F> void handleArray(const O* o, const TObjArray* c, F pmf) {
TObjArrayIter arr(c);
for(TObject* i = arr.Next(); i; i=arr.Next())
(o->*pmf)(i);
}
template <typename O, typename C, typename F> void handleMap(const O* o, const C& c, F pmf) {
for (typename C::const_iterator i = c.begin(); i != c.end(); ++i)
(o->*pmf)((*i).first, (*i).second);
}
template <typename O, typename C, typename F> void handleRMap(const O* o, const C& c, F pmf) {
for (typename C::const_reverse_iterator i = c.rbegin(); i != c.rend(); ++i) {
//cout << "Handle RMAP [ " << (*i).first << " ]" << endl;
handle(o, (*i).second, pmf);
}
template <typename O, typename C, typename F> void handleRMap_(const O* o, const C& c, F pmf) {
for (typename C::const_iterator i = c.begin(); i != c.end(); ++i) {
const auto& cc = (*i).second;
for (const auto& j : cc) {
(o->*pmf)(j);
}
}
}
Markus Frank
committed
Geant4Converter& Geant4Converter::create(DetElement top) {
Geant4GeometryInfo& geo = this->init();
geo.manager = &wrld.detectorDescription().manager();
Markus Frank
committed
collect(top, geo);
Markus Frank
committed
checkOverlaps = false;
// We do not have to handle defines etc.
// All positions and the like are not really named.
// Hence, start creating the G4 objects for materials, solids and log volumes.
Markus Frank
committed
//outputLevel = WARNING;
//debugMaterials = true;
//debugElements = true;
//debugReflections = true;
//debugPlacements = true;
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
handleArray(this, geo.manager->GetListOfGDMLMatrices(), &Geant4Converter::handleMaterialProperties);
handleArray(this, geo.manager->GetListOfOpticalSurfaces(), &Geant4Converter::handleOpticalSurface);
#endif
handle(this, geo.volumes, &Geant4Converter::collectVolume);
handle(this, geo.solids, &Geant4Converter::handleSolid);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "++ Handled %ld solids.", geo.solids.size());
handleRefs(this, geo.vis, &Geant4Converter::handleVis);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "++ Handled %ld visualization attributes.", geo.vis.size());
handleMap(this, geo.limits, &Geant4Converter::handleLimitSet);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "++ Handled %ld limit sets.", geo.limits.size());
handleMap(this, geo.regions, &Geant4Converter::handleRegion);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "++ Handled %ld regions.", geo.regions.size());
handle(this, geo.volumes, &Geant4Converter::handleVolume);
Markus Frank
committed
printout(outputLevel, "Geant4Converter", "++ Handled %ld volumes.", geo.volumes.size());
handleRMap(this, *m_data, &Geant4Converter::handleAssembly);
// Now place all this stuff appropriately
handleRMap(this, *m_data, &Geant4Converter::handlePlacement);
//handleRMap_(this, *m_places, &Geant4Converter::handlePlacement2);
#if ROOT_VERSION_CODE >= ROOT_VERSION(6,17,0)
/// Handle concrete surfaces
handleArray(this, geo.manager->GetListOfSkinSurfaces(), &Geant4Converter::handleSkinSurface);
handleArray(this, geo.manager->GetListOfBorderSurfaces(), &Geant4Converter::handleBorderSurface);
#endif
//==================== Fields
Markus Frank
committed
if ( printSensitives ) {
handleMap(this, geo.sensitives, &Geant4Converter::printSensitive);
}
if ( printPlacements ) {
handleRMap(this, *m_data, &Geant4Converter::printPlacement);
}
geo.setWorld(top.placement().ptr());
geo.valid = true;
printout(INFO, "Geant4Converter", "+++ Successfully converted geometry to Geant4.");
Markus Frank
committed
return *this;