Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-07-08 01:46:14

0001 #include "SimG4Core/Notification/interface/BeginOfRun.h"
0002 #include "SimG4Core/Notification/interface/Observer.h"
0003 #include "SimG4Core/Watcher/interface/SimWatcher.h"
0004 
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "DataFormats/Math/interface/angle_units.h"
0008 
0009 #include "G4Box.hh"
0010 #include "G4Cons.hh"
0011 #include "G4ExtrudedSolid.hh"
0012 #include "G4Polycone.hh"
0013 #include "G4Polyhedra.hh"
0014 #include "G4Trap.hh"
0015 #include "G4Trd.hh"
0016 #include "G4Tubs.hh"
0017 #include "G4LogicalVolumeStore.hh"
0018 #include "G4LogicalVolume.hh"
0019 #include "G4VSolid.hh"
0020 #include "G4NavigationHistory.hh"
0021 #include "G4PhysicalVolumeStore.hh"
0022 #include "G4TransportationManager.hh"
0023 #include "G4VPhysicalVolume.hh"
0024 
0025 #include <algorithm>
0026 #include <map>
0027 #include <set>
0028 #include <sstream>
0029 #include <string>
0030 #include <vector>
0031 
0032 using angle_units::operators::convertRadToDeg;
0033 
0034 class PrintG4Solids : public SimWatcher, public Observer<const BeginOfRun *> {
0035 public:
0036   PrintG4Solids(edm::ParameterSet const &p);
0037   ~PrintG4Solids() override = default;
0038 
0039 private:
0040   void update(const BeginOfRun *run) override;
0041   void dumpSummary(std::ostream &out = G4cout);
0042   G4VPhysicalVolume *getTopPV();
0043 
0044 private:
0045   G4VPhysicalVolume *theTopPV_;
0046 };
0047 
0048 PrintG4Solids::PrintG4Solids(const edm::ParameterSet &p) {
0049   G4cout << "PrintG4Solids:: initialised for printing information about G4VSolids" << G4endl;
0050 }
0051 
0052 void PrintG4Solids::update(const BeginOfRun *run) {
0053   //Now take action
0054   theTopPV_ = getTopPV();
0055 
0056   dumpSummary(G4cout);
0057 }
0058 
0059 void PrintG4Solids::dumpSummary(std::ostream &out) {
0060   //---------- Dump number of objects of each class
0061   out << " @@@@@@@@@@@@@@@@@@ Dumping G4 geometry objects Summary " << G4endl;
0062   if (theTopPV_ == nullptr) {
0063     out << " No volume created " << G4endl;
0064     return;
0065   }
0066   out << " @@@ Geometry built inside world volume: " << theTopPV_->GetName() << G4endl;
0067   // Get number of solids (< # LV if several LV share a solid)
0068   const G4LogicalVolumeStore *lvs = G4LogicalVolumeStore::GetInstance();
0069   std::vector<G4LogicalVolume *>::const_iterator lvcite;
0070   std::set<G4VSolid *> theSolids;
0071   for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
0072     theSolids.insert((*lvcite)->GetSolid());
0073   out << " Number of G4VSolid's: " << theSolids.size() << G4endl;
0074   std::set<G4VSolid *>::const_iterator solid;
0075   for (solid = theSolids.begin(); solid != theSolids.end(); solid++) {
0076     G4String type = (*solid)->GetEntityType();
0077     out << (*solid)->GetName() << ":" << type << " Volume " << (*solid)->GetCubicVolume();
0078     if (type == "G4Box") {
0079       const G4Box *box = static_cast<const G4Box *>(*solid);
0080       out << " dx:dy:dz " << box->GetXHalfLength() << ":" << box->GetYHalfLength() << ":" << box->GetZHalfLength();
0081     } else if (type == "G4Tubs") {
0082       const G4Tubs *tube = static_cast<const G4Tubs *>(*solid);
0083       out << " rin:rout:dz:phistart:dphi " << tube->GetInnerRadius() << ":" << tube->GetOuterRadius() << ":"
0084           << tube->GetZHalfLength() << ":" << convertRadToDeg(tube->GetStartPhiAngle()) << ":"
0085           << convertRadToDeg(tube->GetDeltaPhiAngle());
0086     } else if (type == "G4Cons") {
0087       const G4Cons *cone = static_cast<const G4Cons *>(*solid);
0088       out << " rinminus:routminus:rinplus:routplus:dz:phistart:dphi " << cone->GetInnerRadiusMinusZ() << ":"
0089           << cone->GetOuterRadiusMinusZ() << ":" << cone->GetInnerRadiusPlusZ() << ":" << cone->GetOuterRadiusPlusZ()
0090           << ":" << cone->GetZHalfLength() << ":" << convertRadToDeg(cone->GetStartPhiAngle()) << ":"
0091           << convertRadToDeg(cone->GetDeltaPhiAngle());
0092     } else if (type == "G4Trap") {
0093       const G4Trap *trap = static_cast<const G4Trap *>(*solid);
0094       out << " zhalf:yl1:xl11:xl12:tana1:yl2:xl21:xl22:tana2 " << trap->GetZHalfLength() << ":"
0095           << trap->GetYHalfLength1() << ":" << trap->GetXHalfLength1() << ":" << trap->GetXHalfLength2() << ":"
0096           << trap->GetTanAlpha1() << ":" << trap->GetYHalfLength2() << ":" << trap->GetXHalfLength3() << ":"
0097           << trap->GetXHalfLength4() << ":" << trap->GetTanAlpha2();
0098     } else if (type == "G4Trd") {
0099       const G4Trd *trd = static_cast<const G4Trd *>(*solid);
0100       out << " xl1:xl2:yl1:yl2:zhalf " << trd->GetXHalfLength1() << ":" << trd->GetXHalfLength2() << ":"
0101           << trd->GetYHalfLength1() << ":" << trd->GetYHalfLength2() << ":" << trd->GetZHalfLength();
0102     } else if (type == "G4Polycone") {
0103       const G4Polycone *cone = static_cast<const G4Polycone *>(*solid);
0104       const auto hist = cone->GetOriginalParameters();
0105       int num = hist->Num_z_planes;
0106       out << " angle " << convertRadToDeg(hist->Start_angle) << ":" << convertRadToDeg(hist->Opening_angle) << " with "
0107           << num << " planes:";
0108       for (int k = 0; k < num; ++k)
0109         out << " [" << k << "] " << hist->Z_values[k] << ":" << hist->Rmin[k] << ":" << hist->Rmax[k];
0110     } else if (type == "G4Polyhedra") {
0111       const G4Polyhedra *pgon = static_cast<const G4Polyhedra *>(*solid);
0112       const auto hist = pgon->GetOriginalParameters();
0113       int num = hist->Num_z_planes;
0114       out << " angle " << convertRadToDeg(hist->Start_angle) << ":" << convertRadToDeg(hist->Opening_angle) << " with "
0115           << hist->numSide << " sides and " << num << " planes:";
0116       for (int k = 0; k < num; ++k)
0117         out << " [" << k << "] " << hist->Z_values[k] << ":" << hist->Rmin[k] << ":" << hist->Rmax[k];
0118     } else if (type == "G4ExtrudedSolid") {
0119       const G4ExtrudedSolid *pgon = static_cast<const G4ExtrudedSolid *>(*solid);
0120       int vert = pgon->GetNofVertices();
0121       int numz = pgon->GetNofZSections();
0122       out << " " << vert << " vertices:";
0123       for (int k = 0; k < vert; ++k)
0124         out << " [" << k << "] " << pgon->GetVertex(k);
0125       out << "; and " << numz << " z-sections:";
0126       for (int k = 0; k < numz; ++k) {
0127         const auto &zsec = pgon->GetZSection(k);
0128         out << " [" << k << "] " << zsec.fZ << ":" << zsec.fScale << ":" << zsec.fOffset;
0129       }
0130     }
0131     out << G4endl;
0132   }
0133 }
0134 
0135 G4VPhysicalVolume *PrintG4Solids::getTopPV() {
0136   return G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume();
0137 }
0138 
0139 #include "SimG4Core/Watcher/interface/SimWatcherFactory.h"
0140 #include "FWCore/PluginManager/interface/ModuleDef.h"
0141 
0142 DEFINE_SIMWATCHER(PrintG4Solids);