Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-12-17 09:11:59 

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2020 CERN for the benefit of the Acts project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008 
0009 #include "Acts/Visualization/GeometryView3D.hpp"
0010 
0011 #include "Acts/Detector/DetectorVolume.hpp"
0012 #include "Acts/Detector/Portal.hpp"
0013 #include "Acts/Geometry/BoundarySurfaceT.hpp"
0014 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
0015 #include "Acts/Geometry/Extent.hpp"
0016 #include "Acts/Geometry/GeometryIdentifier.hpp"
0017 #include "Acts/Geometry/Layer.hpp"
0018 #include "Acts/Geometry/Polyhedron.hpp"
0019 #include "Acts/Geometry/TrackingVolume.hpp"
0020 #include "Acts/Geometry/Volume.hpp"
0021 #include "Acts/Surfaces/ConeBounds.hpp"
0022 #include "Acts/Surfaces/ConeSurface.hpp"
0023 #include "Acts/Surfaces/CylinderBounds.hpp"
0024 #include "Acts/Surfaces/CylinderSurface.hpp"
0025 #include "Acts/Surfaces/DiscSurface.hpp"
0026 #include "Acts/Surfaces/RadialBounds.hpp"
0027 #include "Acts/Surfaces/Surface.hpp"
0028 #include "Acts/Surfaces/SurfaceArray.hpp"
0029 #include "Acts/Utilities/BinnedArray.hpp"
0030 #include "Acts/Utilities/BinningType.hpp"
0031 #include "Acts/Utilities/IAxis.hpp"
0032 #include "Acts/Utilities/UnitVectors.hpp"
0033 #include "Acts/Visualization/IVisualization3D.hpp"
0034 
0035 #include <algorithm>
0036 #include <cmath>
0037 #include <memory>
0038 #include <ostream>
0039 #include <utility>
0040 #include <vector>
0041 
0042 #include <limits.h>
0043 #include <unistd.h>
0044 
0045 namespace {
0046 
0047 std::string joinPaths(const std::string& a, const std::string& b) {
0048   if (b.substr(0, 1) == "/" || a.empty()) {
0049     return b;
0050   }
0051 
0052   if (a.substr(a.size() - 1) == "/") {
0053     return a.substr(a.size() - 1) + "/" + b;
0054   }
0055 
0056   return a + "/" + b;
0057 }
0058 
0059 std::string getWorkingDirectory() {
0060   char buffer[PATH_MAX];
0061   return (getcwd(buffer, sizeof(buffer)) != nullptr ? std::string(buffer)
0062                                                     : std::string(""));
0063 }
0064 
0065 }  // namespace
0066 
0067 namespace Acts::Experimental {
0068 ViewConfig s_viewSensitive = ViewConfig({0, 180, 240});
0069 ViewConfig s_viewPassive = ViewConfig({240, 280, 0});
0070 ViewConfig s_viewVolume = ViewConfig({220, 220, 0});
0071 ViewConfig s_viewGrid = ViewConfig({220, 0, 0});
0072 ViewConfig s_viewLine = ViewConfig({0, 0, 220});
0073 }  // namespace Acts::Experimental
0074 
0075 void Acts::GeometryView3D::drawPolyhedron(IVisualization3D& helper,
0076                                           const Polyhedron& polyhedron,
0077                                           const ViewConfig& viewConfig) {
0078   if (viewConfig.visible) {
0079     if (!viewConfig.triangulate) {
0080       helper.faces(polyhedron.vertices, polyhedron.faces, viewConfig.color);
0081     } else {
0082       helper.faces(polyhedron.vertices, polyhedron.triangularMesh,
0083                    viewConfig.color);
0084     }
0085   }
0086 }
0087 
0088 void Acts::GeometryView3D::drawSurface(IVisualization3D& helper,
0089                                        const Surface& surface,
0090                                        const GeometryContext& gctx,
0091                                        const Transform3& transform,
0092                                        const ViewConfig& viewConfig) {
0093   Polyhedron surfaceHedron =
0094       surface.polyhedronRepresentation(gctx, viewConfig.nSegments);
0095   if (!transform.isApprox(Transform3::Identity())) {
0096     surfaceHedron.move(transform);
0097   }
0098   drawPolyhedron(helper, surfaceHedron, viewConfig);
0099 }
0100 
0101 void Acts::GeometryView3D::drawSurfaceArray(
0102     IVisualization3D& helper, const SurfaceArray& surfaceArray,
0103     const GeometryContext& gctx, const Transform3& transform,
0104     const ViewConfig& sensitiveConfig, const ViewConfig& passiveConfig,
0105     const ViewConfig& gridConfig, const std::string& _outputDir) {
0106   std::string outputDir =
0107       _outputDir == "." ? getWorkingDirectory() : _outputDir;
0108   // Draw all the surfaces
0109   Extent arrayExtent;
0110   for (const auto& sf : surfaceArray.surfaces()) {
0111     ViewConfig vConfig = sf->associatedDetectorElement() != nullptr
0112                              ? sensitiveConfig
0113                              : passiveConfig;
0114     drawSurface(helper, *sf, gctx, transform, vConfig);
0115     auto sfExtent = sf->polyhedronRepresentation(gctx, 1).extent();
0116     arrayExtent.extend(sfExtent);
0117   }
0118 
0119   if (!sensitiveConfig.outputName.empty()) {
0120     helper.write(joinPaths(outputDir, sensitiveConfig.outputName));
0121     helper.clear();
0122   }
0123 
0124   double thickness = gridConfig.lineThickness;
0125   // Draw the grid itself
0126   auto binning = surfaceArray.binningValues();
0127   auto axes = surfaceArray.getAxes();
0128   if (!binning.empty() && binning.size() == 2 && axes.size() == 2) {
0129     // Cylinder surface array
0130     if (binning[0] == binPhi && binning[1] == binZ) {
0131       double R = arrayExtent.medium(binR) + gridConfig.offset;
0132       auto phiValues = axes[0]->getBinEdges();
0133       auto zValues = axes[1]->getBinEdges();
0134       ViewConfig gridRadConfig = gridConfig;
0135       gridRadConfig.nSegments = phiValues.size();
0136       // Longitudinal lines
0137       for (auto phi : phiValues) {
0138         double cphi = std::cos(phi);
0139         double sphi = std::sin(phi);
0140         Vector3 p1(R * cphi, R * sphi, axes[1]->getMin());
0141         Vector3 p0(R * cphi, R * sphi, axes[1]->getMax());
0142         drawSegment(helper, transform * p0, transform * p1, gridConfig);
0143       }
0144       CylinderVolumeBounds cvb(R - 0.5 * thickness, R + 0.5 * thickness,
0145                                0.5 * thickness);
0146       auto cvbOrientedSurfaces = cvb.orientedSurfaces();
0147       for (auto z : zValues) {
0148         for (const auto& cvbSf : cvbOrientedSurfaces) {
0149           drawSurface(helper, *cvbSf.surface, gctx,
0150                       Translation3(0., 0., z) * transform, gridRadConfig);
0151         }
0152       }
0153 
0154     } else if (binning[0] == binR && binning[1] == binPhi) {
0155       double z = arrayExtent.medium(binZ) + gridConfig.offset;
0156       auto rValues = axes[0]->getBinEdges();
0157       auto phiValues = axes[1]->getBinEdges();
0158       ViewConfig gridRadConfig = gridConfig;
0159       gridRadConfig.nSegments = phiValues.size();
0160       for (auto r : rValues) {
0161         CylinderVolumeBounds cvb(r - 0.5 * thickness, r + 0.5 * thickness,
0162                                  0.5 * thickness);
0163         auto cvbOrientedSurfaces = cvb.orientedSurfaces();
0164         for (const auto& cvbSf : cvbOrientedSurfaces) {
0165           drawSurface(helper, *cvbSf.surface, gctx,
0166                       Translation3(0., 0., z) * transform, gridRadConfig);
0167         }
0168       }
0169       double rMin = axes[0]->getMin();
0170       double rMax = axes[0]->getMax();
0171       for (auto phi : phiValues) {
0172         double cphi = std::cos(phi);
0173         double sphi = std::sin(phi);
0174         Vector3 p1(rMax * cphi, rMax * sphi, z);
0175         Vector3 p0(rMin * cphi, rMin * sphi, z);
0176         drawSegment(helper, transform * p0, transform * p1, gridConfig);
0177       }
0178     }
0179   }
0180 
0181   if (!gridConfig.outputName.empty()) {
0182     helper.write(joinPaths(outputDir, gridConfig.outputName));
0183     helper.clear();
0184   }
0185 }
0186 
0187 void Acts::GeometryView3D::drawVolume(IVisualization3D& helper,
0188                                       const Volume& volume,
0189                                       const GeometryContext& gctx,
0190                                       const Transform3& transform,
0191                                       const ViewConfig& viewConfig) {
0192   auto bSurfaces = volume.volumeBounds().orientedSurfaces(volume.transform());
0193   for (const auto& bs : bSurfaces) {
0194     drawSurface(helper, *bs.surface, gctx, transform, viewConfig);
0195   }
0196 }
0197 
0198 void Acts::GeometryView3D::drawPortal(IVisualization3D& helper,
0199                                       const Experimental::Portal& portal,
0200                                       const GeometryContext& gctx,
0201                                       const Transform3& transform,
0202                                       const ViewConfig& connected,
0203                                       const ViewConfig& disconnected) {
0204   // color the portal based on if it contains two links(green)
0205   // or one link(red)
0206   auto surface = &(portal.surface());
0207   auto links = &(portal.detectorVolumeUpdaters());
0208   if (links->size() == 2) {
0209     drawSurface(helper, *surface, gctx, transform, connected);
0210   } else {
0211     drawSurface(helper, *surface, gctx, transform, disconnected);
0212   }
0213 }
0214 
0215 void Acts::GeometryView3D::drawDetectorVolume(
0216     IVisualization3D& helper, const Experimental::DetectorVolume& volume,
0217     const GeometryContext& gctx, const Transform3& transform,
0218     const ViewConfig& connected, const ViewConfig& unconnected,
0219     const ViewConfig& viewConfig) {
0220   // draw the surfaces of the mother volume
0221   for (auto surface : volume.surfaces()) {
0222     drawSurface(helper, *surface, gctx, transform, viewConfig);
0223   }
0224 
0225   // draw the envelope first
0226   for (auto portal : volume.portals()) {
0227     drawPortal(helper, *portal, gctx, transform, connected, unconnected);
0228   }
0229 
0230   // recurse if there are subvolumes
0231   for (auto subvolume : volume.volumes()) {
0232     drawDetectorVolume(helper, *subvolume, gctx, transform, connected,
0233                        unconnected, viewConfig);
0234   }
0235 }
0236 
0237 void Acts::GeometryView3D::drawLayer(
0238     IVisualization3D& helper, const Layer& layer, const GeometryContext& gctx,
0239     const ViewConfig& layerConfig, const ViewConfig& sensitiveConfig,
0240     const ViewConfig& gridConfig, const std::string& _outputDir) {
0241   std::string outputDir =
0242       _outputDir == "." ? getWorkingDirectory() : _outputDir;
0243 
0244   if (layerConfig.visible) {
0245     auto layerVolume = layer.representingVolume();
0246     if (layerVolume != nullptr) {
0247       drawVolume(helper, *layerVolume, gctx, Transform3::Identity(),
0248                  layerConfig);
0249     } else {
0250       const auto& layerSurface = layer.surfaceRepresentation();
0251       drawSurface(helper, layerSurface, gctx, Transform3::Identity(),
0252                   layerConfig);
0253     }
0254     if (!layerConfig.outputName.empty()) {
0255       helper.write(joinPaths(outputDir, layerConfig.outputName));
0256       helper.clear();
0257     }
0258   }
0259 
0260   if (sensitiveConfig.visible || gridConfig.visible) {
0261     auto surfaceArray = layer.surfaceArray();
0262     if (surfaceArray != nullptr) {
0263       drawSurfaceArray(helper, *surfaceArray, gctx, Transform3::Identity(),
0264                        sensitiveConfig, layerConfig, gridConfig, outputDir);
0265     }
0266   }
0267 }
0268 
0269 void Acts::GeometryView3D::drawTrackingVolume(
0270     IVisualization3D& helper, const TrackingVolume& tVolume,
0271     const GeometryContext& gctx, const ViewConfig& containerView,
0272     const ViewConfig& volumeView, const ViewConfig& layerView,
0273     const ViewConfig& sensitiveView, const ViewConfig& gridView, bool writeIt,
0274     const std::string& tag, const std::string& _outputDir) {
0275   std::string outputDir =
0276       _outputDir == "." ? getWorkingDirectory() : _outputDir;
0277   if (tVolume.confinedVolumes() != nullptr) {
0278     const auto& subVolumes = tVolume.confinedVolumes()->arrayObjects();
0279     for (const auto& tv : subVolumes) {
0280       drawTrackingVolume(helper, *tv, gctx, containerView, volumeView,
0281                          layerView, sensitiveView, gridView, writeIt, tag,
0282                          outputDir);
0283     }
0284   }
0285 
0286   ViewConfig cConfig = containerView;
0287   ViewConfig vConfig = volumeView;
0288   ViewConfig lConfig = layerView;
0289   ViewConfig sConfig = sensitiveView;
0290   ViewConfig gConfig = gridView;
0291   gConfig.nSegments = 8;
0292 
0293   ViewConfig vcConfig = cConfig;
0294   std::string vname = tVolume.volumeName();
0295   if (writeIt) {
0296     std::vector<std::string> repChar = {"::" /*, "|", " ", "{", "}"*/};
0297     // std::cout << "PRE: " << vname << std::endl;
0298     for (const auto& rchar : repChar) {
0299       while (vname.find(rchar) != std::string::npos) {
0300         vname.replace(vname.find(rchar), rchar.size(), std::string("_"));
0301       }
0302     }
0303     if (tVolume.confinedVolumes() == nullptr) {
0304       vcConfig = vConfig;
0305       vcConfig.outputName = vname + std::string("_boundaries") + tag;
0306     } else {
0307       std::stringstream vs;
0308       vs << "Container";
0309       std::vector<GeometryIdentifier::Value> ids{tVolume.geometryId().volume()};
0310 
0311       for (const auto* current = &tVolume; current->motherVolume() != nullptr;
0312            current = current->motherVolume()) {
0313         ids.push_back(current->motherVolume()->geometryId().volume());
0314       }
0315 
0316       for (std::size_t i = ids.size() - 1; i < ids.size(); --i) {
0317         vs << "_v" << ids[i];
0318       }
0319       vname = vs.str();
0320       vcConfig.outputName = vname + std::string("_boundaries") + tag;
0321     }
0322   }
0323 
0324   auto bSurfaces = tVolume.boundarySurfaces();
0325   for (const auto& bs : bSurfaces) {
0326     drawSurface(helper, bs->surfaceRepresentation(), gctx,
0327                 Transform3::Identity(), vcConfig);
0328   }
0329   if (writeIt) {
0330     std::string outputName = joinPaths(outputDir, vcConfig.outputName);
0331     helper.write(outputName);
0332     helper.clear();
0333   }
0334 
0335   if (tVolume.confinedLayers() != nullptr) {
0336     const auto& layers = tVolume.confinedLayers()->arrayObjects();
0337     std::size_t il = 0;
0338     for (const auto& tl : layers) {
0339       if (writeIt) {
0340         lConfig.outputName =
0341             vname + std::string("_passives_l") + std::to_string(il) + tag;
0342         sConfig.outputName =
0343             vname + std::string("_sensitives_l") + std::to_string(il) + tag;
0344         gConfig.outputName =
0345             vname + std::string("_grids_l") + std::to_string(il) + tag;
0346       }
0347       drawLayer(helper, *tl, gctx, lConfig, sConfig, gConfig, outputDir);
0348       ++il;
0349     }
0350   }
0351 }
0352 
0353 void Acts::GeometryView3D::drawSegmentBase(IVisualization3D& helper,
0354                                            const Vector3& start,
0355                                            const Vector3& end, int arrows,
0356                                            double arrowLength,
0357                                            double arrowWidth,
0358                                            const ViewConfig& viewConfig) {
0359   double thickness = viewConfig.lineThickness;
0360 
0361   // Draw the parameter shaft and cone
0362   auto direction = Vector3(end - start).normalized();
0363   double hlength = 0.5 * Vector3(end - start).norm();
0364 
0365   auto unitVectors = makeCurvilinearUnitVectors(direction);
0366   RotationMatrix3 lrotation;
0367   lrotation.col(0) = unitVectors.first;
0368   lrotation.col(1) = unitVectors.second;
0369   lrotation.col(2) = direction;
0370 
0371   Vector3 lcenter = 0.5 * (start + end);
0372   double alength = (thickness > 0.) ? arrowLength * thickness : 2.;
0373   if (alength > hlength) {
0374     alength = hlength;
0375   }
0376 
0377   if (arrows == 2) {
0378     hlength -= alength;
0379   } else if (arrows != 0) {
0380     hlength -= 0.5 * alength;
0381     lcenter -= Vector3(arrows * 0.5 * alength * direction);
0382   }
0383 
0384   // Line - draw a line
0385   if (thickness > 0.) {
0386     auto ltransform = Transform3::Identity();
0387     ltransform.prerotate(lrotation);
0388     ltransform.pretranslate(lcenter);
0389 
0390     auto lbounds = std::make_shared<CylinderBounds>(thickness, hlength);
0391     auto line = Surface::makeShared<CylinderSurface>(ltransform, lbounds);
0392 
0393     drawSurface(helper, *line, GeometryContext(), Transform3::Identity(),
0394                 viewConfig);
0395   } else {
0396     helper.line(start, end, viewConfig.color);
0397   }
0398 
0399   // Arrowheads - if configured
0400   if (arrows != 0) {
0401     double awith = thickness * arrowWidth;
0402     double alpha = atan2(thickness * arrowWidth, alength);
0403     auto plateBounds = std::make_shared<RadialBounds>(thickness, awith);
0404 
0405     if (arrows > 0) {
0406       auto aetransform = Transform3::Identity();
0407       aetransform.prerotate(lrotation);
0408       aetransform.pretranslate(end);
0409       // Arrow cone
0410       auto coneBounds = std::make_shared<ConeBounds>(alpha, -alength, 0.);
0411       auto cone = Surface::makeShared<ConeSurface>(aetransform, coneBounds);
0412       drawSurface(helper, *cone, GeometryContext(), Transform3::Identity(),
0413                   viewConfig);
0414       // Arrow end plate
0415       auto aptransform = Transform3::Identity();
0416       aptransform.prerotate(lrotation);
0417       aptransform.pretranslate(Vector3(end - alength * direction));
0418 
0419       auto plate = Surface::makeShared<DiscSurface>(aptransform, plateBounds);
0420       drawSurface(helper, *plate, GeometryContext(), Transform3::Identity(),
0421                   viewConfig);
0422     }
0423     if (arrows < 0 || arrows == 2) {
0424       auto astransform = Transform3::Identity();
0425       astransform.prerotate(lrotation);
0426       astransform.pretranslate(start);
0427 
0428       // Arrow cone
0429       auto coneBounds = std::make_shared<ConeBounds>(alpha, 0., alength);
0430       auto cone = Surface::makeShared<ConeSurface>(astransform, coneBounds);
0431       drawSurface(helper, *cone, GeometryContext(), Transform3::Identity(),
0432                   viewConfig);
0433       // Arrow end plate
0434       auto aptransform = Transform3::Identity();
0435       aptransform.prerotate(lrotation);
0436       aptransform.pretranslate(Vector3(start + alength * direction));
0437 
0438       auto plate = Surface::makeShared<DiscSurface>(aptransform, plateBounds);
0439       drawSurface(helper, *plate, GeometryContext(), Transform3::Identity(),
0440                   viewConfig);
0441     }
0442   }
0443 }
0444 
0445 void Acts::GeometryView3D::drawSegment(IVisualization3D& helper,
0446                                        const Vector3& start, const Vector3& end,
0447                                        const ViewConfig& viewConfig) {
0448   drawSegmentBase(helper, start, end, 0, 0., 0., viewConfig);
0449 }
0450 
0451 void Acts::GeometryView3D::drawArrowBackward(
0452     IVisualization3D& helper, const Vector3& start, const Vector3& end,
0453     double arrowLength, double arrowWidth, const ViewConfig& viewConfig) {
0454   drawSegmentBase(helper, start, end, -1, arrowLength, arrowWidth, viewConfig);
0455 }
0456 
0457 void Acts::GeometryView3D::drawArrowForward(
0458     IVisualization3D& helper, const Vector3& start, const Vector3& end,
0459     double arrowLength, double arrowWidth, const ViewConfig& viewConfig) {
0460   drawSegmentBase(helper, start, end, 1, arrowLength, arrowWidth, viewConfig);
0461 }
0462 
0463 void Acts::GeometryView3D::drawArrowsBoth(IVisualization3D& helper,
0464                                           const Vector3& start,
0465                                           const Vector3& end,
0466                                           double arrowLength, double arrowWidth,
0467                                           const ViewConfig& viewConfig) {
0468   drawSegmentBase(helper, start, end, 2, arrowLength, arrowWidth, viewConfig);
0469 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!