src/Visualization/GeometryView3D.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://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/BinningType.hpp"
0030 #include "Acts/Utilities/IAxis.hpp"
0031 #include "Acts/Utilities/UnitVectors.hpp"
0032 #include "Acts/Visualization/IVisualization3D.hpp"
0033
0034 #include <algorithm>
0035 #include <cmath>
0036 #include <filesystem>
0037 #include <memory>
0038 #include <utility>
0039 #include <vector>
0040
0041 void Acts::GeometryView3D::drawPolyhedron(IVisualization3D& helper,
0042                                           const Polyhedron& polyhedron,
0043                                           const ViewConfig& viewConfig) {
0044   polyhedron.visualize(helper, viewConfig);
0045 }
0046
0047 void Acts::GeometryView3D::drawSurface(IVisualization3D& helper,
0048                                        const Surface& surface,
0049                                        const GeometryContext& gctx,
0050                                        const Transform3& /*transform*/,
0051                                        const ViewConfig& viewConfig) {
0052   surface.visualize(helper, gctx, viewConfig);
0053 }
0054
0055 void Acts::GeometryView3D::drawSurfaceArray(
0056     IVisualization3D& helper, const SurfaceArray& surfaceArray,
0057     const GeometryContext& gctx, const Transform3& transform,
0058     const ViewConfig& sensitiveConfig, const ViewConfig& passiveConfig,
0059     const ViewConfig& gridConfig, const std::filesystem::path& outputDir) {
0060   // Draw all the surfaces
0061   Extent arrayExtent;
0062   for (const auto& sf : surfaceArray.surfaces()) {
0063     ViewConfig vConfig = sf->associatedDetectorElement() != nullptr
0064                              ? sensitiveConfig
0065                              : passiveConfig;
0066     drawSurface(helper, *sf, gctx, transform, vConfig);
0067     auto sfExtent = sf->polyhedronRepresentation(gctx, 1).extent();
0068     arrayExtent.extend(sfExtent);
0069   }
0070
0071   if (!sensitiveConfig.outputName.empty()) {
0072     helper.write(outputDir / sensitiveConfig.outputName);
0073     helper.clear();
0074   }
0075
0076   double thickness = gridConfig.lineThickness;
0077   // Draw the grid itself
0078   auto binning = surfaceArray.binningValues();
0079   auto axes = surfaceArray.getAxes();
0080   if (!binning.empty() && binning.size() == 2 && axes.size() == 2) {
0081     // Cylinder surface array
0082     if (binning[0] == AxisDirection::AxisPhi &&
0083         binning[1] == AxisDirection::AxisZ) {
0084       double R = arrayExtent.medium(AxisDirection::AxisR) + gridConfig.offset;
0085       auto phiValues = axes[0]->getBinEdges();
0086       auto zValues = axes[1]->getBinEdges();
0087       ViewConfig gridRadConfig = gridConfig;
0088       // Longitudinal lines
0089       for (auto phi : phiValues) {
0090         double cphi = std::cos(phi);
0091         double sphi = std::sin(phi);
0092         Vector3 p1(R * cphi, R * sphi, axes[1]->getMin());
0093         Vector3 p0(R * cphi, R * sphi, axes[1]->getMax());
0094         drawSegment(helper, transform * p0, transform * p1, gridConfig);
0095       }
0096       CylinderVolumeBounds cvb(R - 0.5 * thickness, R + 0.5 * thickness,
0097                                0.5 * thickness);
0098       auto cvbOrientedSurfaces = cvb.orientedSurfaces();
0099       for (auto z : zValues) {
0100         for (const auto& cvbSf : cvbOrientedSurfaces) {
0101           drawSurface(helper, *cvbSf.surface, gctx,
0102                       Translation3(0., 0., z) * transform, gridRadConfig);
0103         }
0104       }
0105
0106     } else if (binning[0] == AxisDirection::AxisR &&
0107                binning[1] == AxisDirection::AxisPhi) {
0108       double z = arrayExtent.medium(AxisDirection::AxisZ) + gridConfig.offset;
0109       auto rValues = axes[0]->getBinEdges();
0110       auto phiValues = axes[1]->getBinEdges();
0111       ViewConfig gridRadConfig = gridConfig;
0112       gridRadConfig.quarterSegments = phiValues.size();
0113       for (auto r : rValues) {
0114         CylinderVolumeBounds cvb(r - 0.5 * thickness, r + 0.5 * thickness,
0115                                  0.5 * thickness);
0116         auto cvbOrientedSurfaces = cvb.orientedSurfaces();
0117         for (const auto& cvbSf : cvbOrientedSurfaces) {
0118           drawSurface(helper, *cvbSf.surface, gctx,
0119                       Translation3(0., 0., z) * transform, gridRadConfig);
0120         }
0121       }
0122       double rMin = axes[0]->getMin();
0123       double rMax = axes[0]->getMax();
0124       for (auto phi : phiValues) {
0125         double cphi = std::cos(phi);
0126         double sphi = std::sin(phi);
0127         Vector3 p1(rMax * cphi, rMax * sphi, z);
0128         Vector3 p0(rMin * cphi, rMin * sphi, z);
0129         drawSegment(helper, transform * p0, transform * p1, gridConfig);
0130       }
0131     }
0132   }
0133
0134   if (!gridConfig.outputName.empty()) {
0135     helper.write(outputDir / gridConfig.outputName);
0136     helper.clear();
0137   }
0138 }
0139
0140 void Acts::GeometryView3D::drawVolume(IVisualization3D& helper,
0141                                       const Volume& volume,
0142                                       const GeometryContext& gctx,
0143                                       const Transform3& /*transform*/,
0144                                       const ViewConfig& viewConfig) {
0145   volume.visualize(helper, gctx, viewConfig);
0146 }
0147
0148 void Acts::GeometryView3D::drawPortal(IVisualization3D& helper,
0149                                       const Experimental::Portal& portal,
0150                                       const GeometryContext& gctx,
0151                                       const Transform3& transform,
0152                                       const ViewConfig& connected,
0153                                       const ViewConfig& disconnected) {
0154   // color the portal based on if it contains two links(green)
0155   // or one link(red)
0156   auto surface = &(portal.surface());
0157   auto links = &(portal.portalNavigation());
0158   if (links->size() == 2) {
0159     drawSurface(helper, *surface, gctx, transform, connected);
0160   } else {
0161     drawSurface(helper, *surface, gctx, transform, disconnected);
0162   }
0163 }
0164
0165 void Acts::GeometryView3D::drawDetectorVolume(
0166     IVisualization3D& helper, const Experimental::DetectorVolume& volume,
0167     const GeometryContext& gctx, const Transform3& transform,
0168     const ViewConfig& connected, const ViewConfig& unconnected,
0169     const ViewConfig& viewConfig) {
0170   // draw the surfaces of the mother volume
0171   for (auto surface : volume.surfaces()) {
0172     drawSurface(helper, *surface, gctx, transform, viewConfig);
0173   }
0174
0175   // draw the envelope first
0176   for (auto portal : volume.portals()) {
0177     drawPortal(helper, *portal, gctx, transform, connected, unconnected);
0178   }
0179
0180   // recurse if there are subvolumes
0181   for (auto subvolume : volume.volumes()) {
0182     drawDetectorVolume(helper, *subvolume, gctx, transform, connected,
0183                        unconnected, viewConfig);
0184   }
0185 }
0186
0187 void Acts::GeometryView3D::drawLayer(
0188     IVisualization3D& helper, const Layer& layer, const GeometryContext& gctx,
0189     const ViewConfig& layerConfig, const ViewConfig& sensitiveConfig,
0190     const ViewConfig& gridConfig, const std::filesystem::path& outputDir) {
0191   if (layerConfig.visible) {
0192     auto layerVolume = layer.representingVolume();
0193     if (layerVolume != nullptr) {
0194       drawVolume(helper, *layerVolume, gctx, Transform3::Identity(),
0195                  layerConfig);
0196     } else {
0197       const auto& layerSurface = layer.surfaceRepresentation();
0198       drawSurface(helper, layerSurface, gctx, Transform3::Identity(),
0199                   layerConfig);
0200     }
0201     if (!layerConfig.outputName.empty()) {
0202       helper.write(outputDir / layerConfig.outputName);
0203       helper.clear();
0204     }
0205   }
0206
0207   if (sensitiveConfig.visible || gridConfig.visible) {
0208     auto surfaceArray = layer.surfaceArray();
0209     if (surfaceArray != nullptr) {
0210       drawSurfaceArray(helper, *surfaceArray, gctx, Transform3::Identity(),
0211                        sensitiveConfig, layerConfig, gridConfig, outputDir);
0212     }
0213   }
0214 }
0215
0216 void Acts::GeometryView3D::drawTrackingVolume(
0217     IVisualization3D& helper, const TrackingVolume& tVolume,
0218     const GeometryContext& gctx, const ViewConfig& containerView,
0219     const ViewConfig& volumeView, const ViewConfig& layerView,
0220     const ViewConfig& sensitiveView, const ViewConfig& gridView, bool writeIt,
0221     const std::string& tag, const std::filesystem::path& outputDir) {
0222   if (tVolume.confinedVolumes() != nullptr) {
0223     const auto& subVolumes = tVolume.confinedVolumes()->arrayObjects();
0224     for (const auto& tv : subVolumes) {
0225       drawTrackingVolume(helper, *tv, gctx, containerView, volumeView,
0226                          layerView, sensitiveView, gridView, writeIt, tag,
0227                          outputDir);
0228     }
0229   }
0230
0231   ViewConfig cConfig = containerView;
0232   ViewConfig vConfig = volumeView;
0233   ViewConfig lConfig = layerView;
0234   ViewConfig sConfig = sensitiveView;
0235   ViewConfig gConfig = gridView;
0236   gConfig.quarterSegments = 8;
0237
0238   ViewConfig vcConfig = cConfig;
0239   std::string vname = tVolume.volumeName();
0240   if (writeIt) {
0241     std::vector<std::string> repChar = {"::" /*, "|", " ", "{", "}"*/};
0242     for (const auto& rchar : repChar) {
0243       while (vname.find(rchar) != std::string::npos) {
0244         vname.replace(vname.find(rchar), rchar.size(), std::string("_"));
0245       }
0246     }
0247     if (tVolume.confinedVolumes() == nullptr) {
0248       vcConfig = vConfig;
0249       vcConfig.outputName =
0250           std::filesystem::path(vname + std::string("_boundaries") + tag);
0251     } else {
0252       std::vector<GeometryIdentifier::Value> ids{tVolume.geometryId().volume()};
0253
0254       for (const auto* current = &tVolume; current->motherVolume() != nullptr;
0255            current = current->motherVolume()) {
0256         ids.push_back(current->motherVolume()->geometryId().volume());
0257       }
0258
0259       std::ranges::reverse(ids);
0260       vname = "Container";
0261       for (const auto& id : ids) {
0262         vname += "_v" + std::to_string(id);
0263       }
0264
0265       vcConfig.outputName =
0266           std::filesystem::path(vname + std::string("_boundaries") + tag);
0267     }
0268   }
0269
0270   auto bSurfaces = tVolume.boundarySurfaces();
0271   for (const auto& bs : bSurfaces) {
0272     drawSurface(helper, bs->surfaceRepresentation(), gctx,
0273                 Transform3::Identity(), vcConfig);
0274   }
0275   if (writeIt) {
0276     const std::filesystem::path outputName = outputDir / vcConfig.outputName;
0277     helper.write(outputName);
0278     helper.clear();
0279   }
0280
0281   if (tVolume.confinedLayers() != nullptr) {
0282     const auto& layers = tVolume.confinedLayers()->arrayObjects();
0283     std::size_t il = 0;
0284     for (const auto& tl : layers) {
0285       if (writeIt) {
0286         lConfig.outputName = std::filesystem::path(
0287             vname + std::string("_passives_l") + std::to_string(il) + tag);
0288         sConfig.outputName = std::filesystem::path(
0289             vname + std::string("_sensitives_l") + std::to_string(il) + tag);
0290         gConfig.outputName = std::filesystem::path(
0291             vname + std::string("_grids_l") + std::to_string(il) + tag);
0292       }
0293       drawLayer(helper, *tl, gctx, lConfig, sConfig, gConfig, outputDir);
0294       ++il;
0295     }
0296   }
0297 }
0298
0299 void Acts::GeometryView3D::drawSegmentBase(IVisualization3D& helper,
0300                                            const Vector3& start,
0301                                            const Vector3& end, int arrows,
0302                                            double arrowLength,
0303                                            double arrowWidth,
0304                                            const ViewConfig& viewConfig) {
0305   double thickness = viewConfig.lineThickness;
0306
0307   // Draw the parameter shaft and cone
0308   auto direction = Vector3(end - start).normalized();
0309   double hlength = 0.5 * Vector3(end - start).norm();
0310
0311   auto unitVectors = makeCurvilinearUnitVectors(direction);
0312   RotationMatrix3 lrotation;
0313   lrotation.col(0) = unitVectors.first;
0314   lrotation.col(1) = unitVectors.second;
0315   lrotation.col(2) = direction;
0316
0317   Vector3 lcenter = 0.5 * (start + end);
0318   double alength = (thickness > 0.) ? arrowLength * thickness : 2.;
0319   if (alength > hlength) {
0320     alength = hlength;
0321   }
0322
0323   if (arrows == 2) {
0324     hlength -= alength;
0325   } else if (arrows != 0) {
0326     hlength -= 0.5 * alength;
0327     lcenter -= Vector3(arrows * 0.5 * alength * direction);
0328   }
0329
0330   // Line - draw a line
0331   if (thickness > 0.) {
0332     auto ltransform = Transform3::Identity();
0333     ltransform.prerotate(lrotation);
0334     ltransform.pretranslate(lcenter);
0335
0336     auto lbounds = std::make_shared<CylinderBounds>(thickness, hlength);
0337     auto line = Surface::makeShared<CylinderSurface>(ltransform, lbounds);
0338
0339     drawSurface(helper, *line, GeometryContext(), Transform3::Identity(),
0340                 viewConfig);
0341   } else {
0342     helper.line(start, end, viewConfig.color);
0343   }
0344
0345   // Arrowheads - if configured
0346   if (arrows != 0) {
0347     double awith = thickness * arrowWidth;
0348     double alpha = atan2(thickness * arrowWidth, alength);
0349     auto plateBounds = std::make_shared<RadialBounds>(thickness, awith);
0350
0351     if (arrows > 0) {
0352       auto aetransform = Transform3::Identity();
0353       aetransform.prerotate(lrotation);
0354       aetransform.pretranslate(end);
0355       // Arrow cone
0356       auto coneBounds = std::make_shared<ConeBounds>(alpha, -alength, 0.);
0357       auto cone = Surface::makeShared<ConeSurface>(aetransform, coneBounds);
0358       drawSurface(helper, *cone, GeometryContext(), Transform3::Identity(),
0359                   viewConfig);
0360       // Arrow end plate
0361       auto aptransform = Transform3::Identity();
0362       aptransform.prerotate(lrotation);
0363       aptransform.pretranslate(Vector3(end - alength * direction));
0364
0365       auto plate = Surface::makeShared<DiscSurface>(aptransform, plateBounds);
0366       drawSurface(helper, *plate, GeometryContext(), Transform3::Identity(),
0367                   viewConfig);
0368     }
0369     if (arrows < 0 || arrows == 2) {
0370       auto astransform = Transform3::Identity();
0371       astransform.prerotate(lrotation);
0372       astransform.pretranslate(start);
0373
0374       // Arrow cone
0375       auto coneBounds = std::make_shared<ConeBounds>(alpha, 0., alength);
0376       auto cone = Surface::makeShared<ConeSurface>(astransform, coneBounds);
0377       drawSurface(helper, *cone, GeometryContext(), Transform3::Identity(),
0378                   viewConfig);
0379       // Arrow end plate
0380       auto aptransform = Transform3::Identity();
0381       aptransform.prerotate(lrotation);
0382       aptransform.pretranslate(Vector3(start + alength * direction));
0383
0384       auto plate = Surface::makeShared<DiscSurface>(aptransform, plateBounds);
0385       drawSurface(helper, *plate, GeometryContext(), Transform3::Identity(),
0386                   viewConfig);
0387     }
0388   }
0389 }
0390
0391 void Acts::GeometryView3D::drawSegment(IVisualization3D& helper,
0392                                        const Vector3& start, const Vector3& end,
0393                                        const ViewConfig& viewConfig) {
0394   drawSegmentBase(helper, start, end, 0, 0., 0., viewConfig);
0395 }
0396
0397 void Acts::GeometryView3D::drawArrowBackward(
0398     IVisualization3D& helper, const Vector3& start, const Vector3& end,
0399     double arrowLength, double arrowWidth, const ViewConfig& viewConfig) {
0400   drawSegmentBase(helper, start, end, -1, arrowLength, arrowWidth, viewConfig);
0401 }
0402
0403 void Acts::GeometryView3D::drawArrowForward(
0404     IVisualization3D& helper, const Vector3& start, const Vector3& end,
0405     double arrowLength, double arrowWidth, const ViewConfig& viewConfig) {
0406   drawSegmentBase(helper, start, end, 1, arrowLength, arrowWidth, viewConfig);
0407 }
0408
0409 void Acts::GeometryView3D::drawArrowsBoth(IVisualization3D& helper,
0410                                           const Vector3& start,
0411                                           const Vector3& end,
0412                                           double arrowLength, double arrowWidth,
0413                                           const ViewConfig& viewConfig) {
0414   drawSegmentBase(helper, start, end, 2, arrowLength, arrowWidth, viewConfig);
0415 }