include/root/TGeoParallelWorld.h

0001 /*************************************************************************
0002  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
0003  * All rights reserved.                                                  *
0004  *                                                                       *
0005  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0006  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0007  *************************************************************************/
0008
0009 // Authors: Andrei Gheata   30/06/14
0010 //          Sandro Wenzel   01/09/24
0011
0012 #ifndef ROOT_TGeoParallelWorld
0013 #define ROOT_TGeoParallelWorld
0014
0015 #include "TNamed.h"
0016 #include "TGeoVoxelGrid.h"
0017
0018 // forward declarations
0019 class TGeoManager;
0020 class TGeoPhysicalNode;
0021 class TGeoVolume;
0022
0023 class TGeoParallelWorld : public TNamed {
0024
0025 public:
0026    // internal enum letting choose between
0027    // VoxelFinder or BVH-based algorithms
0028    enum class AccelerationMode { kVoxelFinder, kBVH };
0029
0030    // a structure for safety evaluation (caching) purpose
0031    // to be stored per 3D grid voxel
0032    struct SafetyVoxelInfo {
0033       float min_safety{-1.f}; // the minimum safety from the mid-point of this voxel to any leaf bounding box
0034       int idx_start{-1}; // the index into an external storage, where candidate bounding boxes to search for this voxel
0035                          // are stored (if -1) --> VoxelInfo not yet initialized
0036       unsigned int num_candidates{0}; // the number of candidate bounding boxes to search
0037       bool isInitialized() const { return idx_start >= 0; }
0038    };
0039
0040 protected:
0041    TGeoManager *fGeoManager;     // base geometry
0042    TObjArray *fPaths;            // array of paths
0043    Bool_t fUseOverlaps;          // Activated if user defined overlapping candidates
0044    Bool_t fIsClosed;             //! Closed flag
0045    TGeoVolume *fVolume;          //! helper volume
0046    TGeoPhysicalNode *fLastState; //! Last PN touched
0047    TObjArray *fPhysical;         //! array of physical nodes
0048
0049    void *fBVH = nullptr; //! BVH helper structure for safety and navigation
0050    TGeoVoxelGrid<SafetyVoxelInfo> *fSafetyVoxelCache =
0051       nullptr;                                        //! A regular 3D cache layer for fast point-based safety lookups
0052    std::vector<unsigned int> fSafetyCandidateStore{}; //! stores bounding boxes serving a quick safety candidates (to be
0053                                                       //! used with the VoxelGrid and SafetyVoxelInfo)
0054    void *fBoundingBoxes = nullptr;                    //! to keep the vector of primitive axis aligned bounding boxes
0055    AccelerationMode fAccMode = AccelerationMode::kVoxelFinder; //! switch between different algorithm implementations
0056
0057    TGeoParallelWorld(const TGeoParallelWorld &) = delete;
0058    TGeoParallelWorld &operator=(const TGeoParallelWorld &) = delete;
0059
0060 public:
0061    // constructors
0062    TGeoParallelWorld()
0063       : TNamed(),
0064         fGeoManager(nullptr),
0065         fPaths(nullptr),
0066         fUseOverlaps(kFALSE),
0067         fIsClosed(kFALSE),
0068         fVolume(nullptr),
0069         fLastState(nullptr),
0070         fPhysical(nullptr)
0071    {
0072    }
0073    TGeoParallelWorld(const char *name, TGeoManager *mgr);
0074
0075    // destructor
0076    ~TGeoParallelWorld() override;
0077    // API for adding components nodes
0078    void AddNode(const char *path);
0079    // Activate/deactivate  overlap usage
0080    void SetUseOverlaps(Bool_t flag) { fUseOverlaps = flag; }
0081    Bool_t IsUsingOverlaps() const { return fUseOverlaps; }
0082    void ResetOverlaps() const;
0083    // Adding overlap candidates can highly improve performance.
0084    void AddOverlap(TGeoVolume *vol, Bool_t activate = kTRUE);
0085    void AddOverlap(const char *volname, Bool_t activate = kTRUE);
0086    // The normal PW mode (without declaring overlaps) does detect them
0087    Int_t PrintDetectedOverlaps() const;
0088
0089    // Closing a parallel geometry is mandatory
0090    Bool_t CloseGeometry();
0091    // Refresh structures in case of re-alignment
0092    void RefreshPhysicalNodes();
0093
0094    // ability to choose algorithm implementation; should be called before CloseGeometry
0095    void SetAccelerationMode(AccelerationMode const &mode) { fAccMode = mode; }
0096    AccelerationMode const &GetAccelerationMode() const { return fAccMode; }
0097
0098    // BVH related functions for building, printing, checking
0099    void BuildBVH();
0100    void PrintBVH() const;
0101    bool CheckBVH(void *, size_t) const;
0102
0103    // --- main navigation interfaces ----
0104
0105    // FindNode
0106    TGeoPhysicalNode *FindNode(Double_t point[3])
0107    {
0108       switch (fAccMode) {
0109       case AccelerationMode::kVoxelFinder: return FindNodeOrig(point);
0110       case AccelerationMode::kBVH:
0111          return FindNodeBVH(point);
0112          // case AccelerationMode::kLoop : return FindNodeLoop(point);
0113       }
0114       return nullptr;
0115    }
0116
0117    // main interface for Safety
0118    Double_t Safety(Double_t point[3], Double_t safmax = 1.E30)
0119    {
0120       switch (fAccMode) {
0121       case AccelerationMode::kVoxelFinder: return SafetyOrig(point, safmax);
0122       case AccelerationMode::kBVH:
0123          return VoxelSafety(point, safmax);
0124          // case AccelerationMode::kLoop : return SafetyLoop(point, safmax);
0125       }
0126       return 0;
0127    }
0128
0129    // main interface for FindNextBoundary
0130    TGeoPhysicalNode *FindNextBoundary(Double_t point[3], Double_t dir[3], Double_t &step, Double_t stepmax = 1.E30)
0131    {
0132       switch (fAccMode) {
0133       case AccelerationMode::kVoxelFinder: return FindNextBoundaryOrig(point, dir, step, stepmax);
0134       case AccelerationMode::kBVH:
0135          return FindNextBoundaryBVH(point, dir, step, stepmax);
0136          // case AccelerationMode::kLoop : return FindNextBoundaryLoop(point, dir, step, stepmax);
0137       }
0138       return nullptr;
0139    }
0140
0141    // Getters
0142    TGeoManager *GetGeometry() const { return fGeoManager; }
0143    Bool_t IsClosed() const { return fIsClosed; }
0144    TGeoVolume *GetVolume() const { return fVolume; }
0145
0146    // Utilities
0147    void CheckOverlaps(Double_t ovlp = 0.001); // default 10 microns
0148    void Draw(Option_t *option) override;
0149
0150 private:
0151    // various implementations for FindNextBoundary
0152    TGeoPhysicalNode *FindNextBoundaryLoop(Double_t point[3], Double_t dir[3], Double_t &step, Double_t stepmax = 1.E30);
0153    TGeoPhysicalNode *FindNextBoundaryOrig(Double_t point[3], Double_t dir[3], Double_t &step, Double_t stepmax = 1.E30);
0154    TGeoPhysicalNode *FindNextBoundaryBVH(Double_t point[3], Double_t dir[3], Double_t &step, Double_t stepmax = 1.E30);
0155
0156    // various implementations for FindNode
0157    TGeoPhysicalNode *FindNodeLoop(Double_t point[3]);
0158    TGeoPhysicalNode *FindNodeOrig(Double_t point[3]);
0159    TGeoPhysicalNode *FindNodeBVH(Double_t point[3]);
0160
0161    // various implementations for Safety
0162    Double_t SafetyLoop(Double_t point[3], Double_t safmax = 1.E30);
0163    Double_t SafetyBVH(Double_t point[3], Double_t safmax = 1.E30);
0164    Double_t SafetyOrig(Double_t point[3], Double_t safmax = 1.E30);
0165    Double_t VoxelSafety(Double_t point[3], Double_t safmax = 1.E30);
0166
0167    // helper functions related to local safety caching (3D voxel grid)
0168
0169    // Given a 3D point, returns the
0170    // a) the min radius R such that there is at least one leaf bounding box fully enclosed
0171    //    in the sphere of radius R around point
0172    // b) the set of leaf bounding boxes who partially lie within radius + margin
0173
0174    // ... using BVH
0175    std::pair<double, double>
0176    GetBVHSafetyCandidates(double point[3], std::vector<int> &candidates, double margin = 0.) const;
0177    // .... same with a simpler, slower algorithm
0178    std::pair<double, double>
0179    GetLoopSafetyCandidates(double point[3], std::vector<int> &candidates, double margin = 0.) const;
0180
0181    void InitSafetyVoxel(TGeoVoxelGridIndex const &);
0182    void TestVoxelGrid(); // a debug method to play with the voxel grid
0183
0184    ClassDefOverride(TGeoParallelWorld, 3) // parallel world base class
0185 };
0186
0187 #endif