Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-05-20 02:50:32

0001 #include "FWCore/Utilities/interface/Exception.h"
0002 
0003 #include "SimG4Core/Geometry/interface/DDG4Builder.h"
0004 #include "SimG4Core/Geometry/interface/DDG4SensitiveConverter.h"
0005 #include "SimG4Core/Geometry/interface/DDG4SolidConverter.h"
0006 #include "SimG4Core/Geometry/interface/SensitiveDetectorCatalog.h"
0007 
0008 #include "DetectorDescription/Core/interface/DDCompactView.h"
0009 #include "DetectorDescription/Core/interface/DDSpecifics.h"
0010 
0011 #include "G4LogicalVolume.hh"
0012 #include "G4Material.hh"
0013 #include "G4PVPlacement.hh"
0014 #include "G4ReflectionFactory.hh"
0015 #include "G4VPhysicalVolume.hh"
0016 #include "G4VSolid.hh"
0017 
0018 #include "G4SystemOfUnits.hh"
0019 #include "G4UnitsTable.hh"
0020 
0021 #include <sstream>
0022 
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 
0025 DDG4Builder::DDG4Builder(const DDCompactView *cpv, G4LogicalVolumeToDDLogicalPartMap &lvmap, bool check)
0026     : solidConverter_(new DDG4SolidConverter), compactView_(cpv), map_(lvmap), check_(check) {
0027   theVectorOfDDG4Dispatchables_ = new DDG4DispContainer();
0028 }
0029 
0030 DDG4Builder::~DDG4Builder() { delete solidConverter_; }
0031 
0032 G4LogicalVolume *DDG4Builder::convertLV(const DDLogicalPart &part) {
0033   edm::LogVerbatim("SimG4CoreGeometry") << "DDG4Builder::convertLV(): DDLogicalPart = " << part;
0034   G4LogicalVolume *result = logs_[part];
0035   if (!result) {
0036     G4VSolid *g4s = convertSolid(part.solid());
0037     G4Material *g4m = convertMaterial(part.material());
0038     result = new G4LogicalVolume(g4s, g4m, part.name().name());
0039     map_.insert(result, part);
0040     DDG4Dispatchable *disp = new DDG4Dispatchable(&part, result);
0041     theVectorOfDDG4Dispatchables_->push_back(disp);
0042     edm::LogVerbatim("SimG4CoreGeometry")
0043         << "DDG4Builder::convertLV(): new G4LogicalVolume " << part.name().name()
0044         << "\nDDG4Builder: newEvent: dd=" << part.ddname() << " g4=" << result->GetName();
0045     logs_[part] = result;  // DDD -> GEANT4
0046   }
0047   return result;
0048 }
0049 
0050 G4VSolid *DDG4Builder::convertSolid(const DDSolid &solid) {
0051   G4VSolid *result = sols_[solid];
0052   if (!result) {
0053     result = solidConverter_->convert(solid);
0054     sols_[solid] = result;
0055   }
0056   return result;
0057 }
0058 
0059 G4Material *DDG4Builder::convertMaterial(const DDMaterial &material) {
0060   edm::LogVerbatim("SimG4CoreGeometry") << "DDDetConstr::ConvertMaterial: material=" << material;
0061   G4Material *result = nullptr;
0062   if (material) {
0063     // only if it's a valid DDD-material
0064     if ((result = mats_[material])) {
0065       edm::LogVerbatim("SimG4CoreGeometry") << "  is already converted";
0066       return result;
0067     }
0068   } else {
0069     // only if it's NOT a valid DDD-material
0070     throw cms::Exception("SimG4CoreGeometry",
0071                          " material is not valid from the Detector Description: " + material.toString());
0072   }
0073   int c = 0;
0074   if ((c = material.noOfConstituents())) {
0075     // it's a composite material
0076     edm::LogVerbatim("SimG4CoreGeometry")
0077         << "  creating a G4-composite material. c=" << c << " d=" << material.density() / CLHEP::g * CLHEP::mole;
0078     result = new G4Material(material.name().name(), material.density(), c);
0079     for (int i = 0; i < c; ++i) {
0080       // recursive building of constituents
0081       edm::LogVerbatim("SimG4CoreGeometry")
0082           << "  adding the composite=" << material.name() << " fm=" << material.constituent(i).second;
0083       result->AddMaterial(convertMaterial(material.constituent(i).first),
0084                           material.constituent(i).second);  // fractionmass
0085     }
0086   } else {
0087     // it's an elementary material
0088     edm::LogVerbatim("SimG4CoreGeometry") << "  building an elementary material"
0089                                           << " z=" << material.z() << " a=" << material.a() / CLHEP::g * CLHEP::mole
0090                                           << " d=" << material.density() / CLHEP::g * CLHEP::cm3;
0091     result = new G4Material(material.name().name(), material.z(), material.a(), material.density());
0092   }
0093   mats_[material] = result;
0094   return result;
0095 }
0096 
0097 G4LogicalVolume *DDG4Builder::BuildGeometry(SensitiveDetectorCatalog &catalog) {
0098   G4ReflectionFactory *refFact = G4ReflectionFactory::Instance();
0099   refFact->SetScalePrecision(100. * refFact->GetScalePrecision());
0100 
0101   using Graph = DDCompactView::Graph;
0102   const auto &gra = compactView_->graph();
0103   using adjl_iterator = Graph::const_adj_iterator;
0104   adjl_iterator git = gra.begin();
0105   adjl_iterator gend = gra.end();
0106 
0107   Graph::index_type i = 0;
0108   for (; git != gend; ++git) {
0109     const DDLogicalPart &ddLP = gra.nodeData(git);
0110     if (!(ddLP.isDefined().second)) {
0111       throw cms::Exception("SimG4CoreGeometry",
0112                            " DDG4Builder::BuildGeometry() has encountered an "
0113                            "undefined DDLogicalPart named " +
0114                                ddLP.toString());
0115     }
0116     G4LogicalVolume *g4LV = convertLV(ddLP);
0117     ++i;
0118     if (!git->empty()) {
0119       // ask for children of ddLP
0120       Graph::edge_list::const_iterator cit = git->begin();
0121       Graph::edge_list::const_iterator cend = git->end();
0122       for (; cit != cend; ++cit) {
0123         // fetch specific data
0124         const DDLogicalPart &ddcurLP = gra.nodeData(cit->first);
0125         if (!ddcurLP.isDefined().second) {
0126           std::string err = " DDG4Builder::BuildGeometry() in processing \"children\" has ";
0127           err += "encountered an undefined DDLogicalPart named " + ddcurLP.toString() + " is a child of " +
0128                  ddLP.toString();
0129           throw cms::Exception("SimG4CoreGeometry", err);
0130         }
0131         int offset = getInt("CopyNoOffset", ddcurLP);
0132         int tag = getInt("CopyNoTag", ddcurLP);
0133         DDRotationMatrix rm(gra.edgeData(cit->second)->rot());
0134         DD3Vector x, y, z;
0135         rm.GetComponents(x, y, z);
0136         if ((x.Cross(y)).Dot(z) < 0)
0137           edm::LogVerbatim("SimG4CoreGeometry")
0138               << "DDG4Builder: Reflection: " << gra.edgeData(cit->second)->ddrot()
0139               << ">>Placement d=" << gra.nodeData(cit->first).ddname() << " m=" << ddLP.ddname()
0140               << " cp=" << gra.edgeData(cit->second)->copyno() << " r=" << gra.edgeData(cit->second)->ddrot().ddname();
0141         G4ThreeVector tempTran(gra.edgeData(cit->second)->trans().X(),
0142                                gra.edgeData(cit->second)->trans().Y(),
0143                                gra.edgeData(cit->second)->trans().Z());
0144         G4Translate3D transl = tempTran;
0145         CLHEP::HepRep3x3 temp(x.X(), x.Y(), x.Z(), y.X(), y.Y(), y.Z(), z.X(), z.Y(), z.Z());  // matrix
0146         CLHEP::HepRotation hr(temp);
0147         edm::LogVerbatim("SimG4CoreGeometry")
0148             << "Position " << gra.nodeData(cit->first).name().name() << ":"
0149             << gra.edgeData(cit->second)->copyno() + offset + tag << " in " << g4LV->GetName() << " at " << tempTran
0150             << " with rotation matrix (" << x.X() << ", " << x.Y() << ", " << x.Z() << ", " << y.X() << ", " << y.Y()
0151             << ", " << y.Z() << ", " << z.X() << ", " << z.Y() << ", " << z.Z() << ")";
0152 
0153         // G3 convention of defining rot-matrices ...
0154         G4Transform3D trfrm = transl * G4Rotate3D(hr.inverse());  //.inverse();
0155 
0156         refFact->Place(trfrm,  // transformation containing a possible reflection
0157                        gra.nodeData(cit->first).name().name(),
0158                        convertLV(gra.nodeData(cit->first)),                 // daugther
0159                        g4LV,                                                // mother
0160                        false,                                               // 'ONLY'
0161                        gra.edgeData(cit->second)->copyno() + offset + tag,  // copy number
0162                        check_);
0163       }  // iterate over children
0164     }    // if (children)
0165   }      // iterate over graph nodes
0166 
0167   // Looking for in the G4ReflectionFactory secretly created reflected
0168   // G4LogicalVolumes
0169   std::map<DDLogicalPart, G4LogicalVolume *>::const_iterator ddg4_it = logs_.begin();
0170   for (; ddg4_it != logs_.end(); ++ddg4_it) {
0171     G4LogicalVolume *reflLogicalVolume = refFact->GetReflectedLV(ddg4_it->second);
0172     if (reflLogicalVolume) {
0173       DDLogicalPart ddlv = ddg4_it->first;
0174       map_.insert(reflLogicalVolume, ddlv);
0175       DDG4Dispatchable *disp = new DDG4Dispatchable(&(ddg4_it->first), reflLogicalVolume);
0176       theVectorOfDDG4Dispatchables_->push_back(disp);
0177       edm::LogVerbatim("SimG4CoreGeometry")
0178           << "DDG4Builder: dd=" << ddlv.ddname() << " g4=" << reflLogicalVolume->GetName();
0179     }
0180   }
0181 
0182   G4LogicalVolume *world = logs_[compactView_->root()];
0183 
0184   //
0185   //  needed for building sensitive detectors
0186   //
0187   DDG4SensitiveConverter conv;
0188   conv.upDate(*theVectorOfDDG4Dispatchables_, catalog);
0189 
0190   return world;
0191 }
0192 
0193 int DDG4Builder::getInt(const std::string &ss, const DDLogicalPart &part) {
0194   DDValue val(ss);
0195   std::vector<const DDsvalues_type *> result = part.specifics();
0196   bool foundIt = false;
0197   for (auto stype : result) {
0198     foundIt = DDfetch(stype, val);
0199     if (foundIt)
0200       break;
0201   }
0202   if (foundIt) {
0203     std::vector<double> temp = val.doubles();
0204     if (temp.size() != 1) {
0205       throw cms::Exception("SimG4CoreGeometry",
0206                            " DDG4Builder::getInt() Problem with Region tags - "
0207                            "one and only one allowed: " +
0208                                ss);
0209     }
0210     return int(temp[0]);
0211   } else
0212     return 0;
0213 }
0214 
0215 double DDG4Builder::getDouble(const std::string &ss, const DDLogicalPart &part) {
0216   DDValue val(ss);
0217   std::vector<const DDsvalues_type *> result = part.specifics();
0218   bool foundIt = false;
0219   for (auto stype : result) {
0220     foundIt = DDfetch(stype, val);
0221     if (foundIt)
0222       break;
0223   }
0224   if (foundIt) {
0225     std::vector<std::string> temp = val.strings();
0226     if (temp.size() != 1) {
0227       throw cms::Exception("SimG4CoreGeometry",
0228                            " DDG4Builder::getDouble() Problem with Region tags "
0229                            "- one and only one allowed: " +
0230                                ss);
0231     }
0232     double v;
0233     std::string unit;
0234     std::istringstream is(temp[0]);
0235     is >> v >> unit;
0236     v = v * G4UnitDefinition::GetValueOf(unit.substr(1, unit.size()));
0237     return v;
0238   } else
0239     return 0;
0240 }