diff --git a/include/xdg/gprt/ray.h b/include/xdg/gprt/ray.h new file mode 100644 index 00000000..023b1338 --- /dev/null +++ b/include/xdg/gprt/ray.h @@ -0,0 +1,42 @@ +#ifndef _XDG_GPRT_RAY_H +#define _XDG_GPRT_RAY_H + +#include "gprt.h" +#include "../shared_enums.h" + +/* + * Double-precision ray and hit structures used by the GPRT backend. + * + * These types are not inherently GPRT-specific, but we keep them here for now + * since GPRT is the only GPU backend. If another GPU backend is added, these + * can be reused. Unifying them with the CPU/Embree types is possible, but may + * not be worth the added complexity at this stage. + */ + +namespace xdg { + +struct dblRay +{ + double3 origin; + double3 direction; + int volume_mesh_id; // MeshID of the volume this ray will be traced against + uint enabled; // Flag to indicate if the ray is active + // TODO - Implement exclude primtives functionality. Right now these are essentially just stubs. + int32_t* exclude_primitives; // Optional for excluding primitives + int32_t exclude_count; // Number of excluded primitives +}; + + +// TODO - Should we define separate hit structs for PIV and ray-fire or do we think its better to keep them together? +struct dblHit +{ + double distance; + int surf_id; + int primitive_id; + PointInVolume piv; // Point in volume check result (0 for outside, 1 for inside) +}; + +} + + +#endif diff --git a/include/xdg/gprt/ray_tracer.h b/include/xdg/gprt/ray_tracer.h index 8d24d107..37e47898 100644 --- a/include/xdg/gprt/ray_tracer.h +++ b/include/xdg/gprt/ray_tracer.h @@ -7,12 +7,9 @@ #include "xdg/constants.h" #include "xdg/mesh_manager_interface.h" -#include "xdg/primitive_ref.h" -#include "xdg/geometry_data.h" #include "xdg/ray_tracing_interface.h" -#include "xdg/ray.h" #include "xdg/error.h" -#include "gprt/gprt.h" + #include "shared_structs.h" extern GPRTProgram dbl_deviceCode; @@ -26,17 +23,16 @@ enum class RayGenType { }; struct gprtRayHit { - size_t capacity = 1; // Max number of rays allocated - size_t size = 0; // Current number of active rays + DeviceRayHitBuffers view; // external facing POD for rayhit buffers + size_t size = 0; // Current number of active rays GPRTBufferOf ray = nullptr; GPRTBufferOf hit = nullptr; - dblRay* devRayAddr = nullptr; - dblHit* devHitAddr = nullptr; - bool is_valid() const { return capacity > 0 && ray && hit && devRayAddr && devHitAddr; } + bool is_valid() const { + return view.capacity > 0 && ray && hit && view.rayDevPtr && view.hitDevPtr; + } }; - class GPRTRayTracer : public RayTracer { public: GPRTRayTracer(); @@ -90,8 +86,17 @@ class GPRTRayTracer : public RayTracer { HitOrientation orientation = HitOrientation::EXITING, std::vector* const exclude_primitives = nullptr) override; + void ray_fire_prepared(const size_t num_rays, + const double dist_limit = INFTY, + HitOrientation orientation = HitOrientation::EXITING) override; + + void point_in_volume_prepared(const size_t num_rays) override; + std::pair closest(TreeID scene, - const Position& origin) override {}; + const Position& origin) override { + fatal_error("Closest queries are not currently supported with GPRT ray tracer"); + return {INFTY, ID_NONE}; + }; bool occluded(TreeID scene, const Position& origin, @@ -100,9 +105,31 @@ class GPRTRayTracer : public RayTracer { fatal_error("Occlusion queries are not currently supported with GPRT ray tracer"); return false; } - + + // Check to see if buffers large enough and resize if not + void check_rayhit_buffer_capacity(const size_t N) override; + + // Method to expose device ray and hit buffers for external population + DeviceRayHitBuffers get_device_rayhit_buffers(const size_t N) override; + + /** + * @brief Allocate device buffers and invoke a callback to populate them + * + * This method enables downstream applications to populate ray buffers using + * any compute API (GPRT, CUDA, HIP, etc.) without XDG needing to know the details. + */ + void populate_rays_external(size_t numRays, + const RayPopulationCallback& callback) override; + + void transfer_hits_buffer_to_host(const size_t num_rays, + std::vector& hits); + + GPRTContext context() + { + return context_; + } + private: - void check_ray_buffer_capacity(size_t N); // GPRT objects GPRTContext context_; @@ -111,9 +138,8 @@ class GPRTRayTracer : public RayTracer { GPRTAccel world_; GPRTBuildParams buildParams_; //> rayGenPrograms_; - GPRTMissOf missProgram_; GPRTComputeOf aabbPopulationProgram_; // surface_volume_tree_to_accel_map; // Map from XDG::TreeID to GPRTAccel for volume TLAS - std::vector blas_handles_; // Store BLAS handles so that they can be explicitly referenced in destructor + std::unordered_map surface_tree_to_volume_map_; + std::vector tlas_handles_; // Host side storage of TLAS device addresses + GPRTBufferOf tlas_handle_buffer_; // Device buffer for TLAS addresses + std::vector meshid_to_sense_; // Host-side MeshID -> sense map + GPRTBufferOf meshid_to_sense_buffer_ {nullptr}; // Device buffer for MeshID -> sense map + bool initialized_ {false}; // flag to indicate if init() has been called + + void update_tlas_table_(); + void update_meshid_to_sense_(); + + // Internal GPRT helper method to upload data to device buffers, creating or resizing as needed + template + void upload_device_buffer_(GPRTBufferOf& buf, const std::vector& host_data) + { + if (host_data.empty()) return; + + if (!buf) { + buf = gprtDeviceBufferCreate(context_, host_data.size(), host_data.data()); + return; + } + + gprtBufferResize(context_, buf, host_data.size(), false); + gprtBufferMap(buf); + std::copy(host_data.begin(), host_data.end(), gprtBufferGetHostPointer(buf)); + gprtBufferUnmap(buf); + } // Global Tree IDs GPRTAccel global_surface_accel_ {nullptr}; GPRTAccel global_element_accel_ {nullptr}; - }; - } // namespace xdg - #endif // include guard \ No newline at end of file diff --git a/include/xdg/gprt/shared_structs.h b/include/xdg/gprt/shared_structs.h index 3855775a..1cff1751 100644 --- a/include/xdg/gprt/shared_structs.h +++ b/include/xdg/gprt/shared_structs.h @@ -1,32 +1,16 @@ +#ifndef XDG_GPRT_SHARED_STRUCTS_H +#define XDG_GPRT_SHARED_STRUCTS_H + #include "gprt.h" #include "../shared_enums.h" +#include "ray.h" struct GPRTPrimitiveRef { int id; // ID of the primitive - int sense; -}; - -struct dblRay -{ - double3 origin; - double3 direction; - double tMin; // Minimum distance for ray intersection - double tMax; // Maximum distance for ray intersection - int32_t* exclude_primitives; // Optional for excluding primitives - int32_t exclude_count; // Number of excluded primitives - xdg::HitOrientation hitOrientation; - int volume_tree; // TreeID of the volume being queried - SurfaceAccelerationStructure volume_accel; // The volume accel + // TODO - What else do we need here? Perhaps a flag for exclude prims? }; -struct dblHit -{ - double distance; - int surf_id; - int primitive_id; - xdg::PointInVolume piv; // Point in volume check result (0 for outside, 1 for inside) -}; /* variables for double precision triangle mesh geometry */ struct DPTriangleGeomData { @@ -35,26 +19,24 @@ struct DPTriangleGeomData { uint3 *index; // index buffer double3 *normals; // normals buffer int surf_id; - int2 vols; - int forward_vol; - int reverse_vol; - dblRay *ray; // double precision rays - xdg::HitOrientation hitOrientation; - int forward_tree; // TreeID of the forward volume - int reverse_tree; // TreeID of the reverse volume + int* meshid_to_sense; // MeshID -> sense (+1 forward, -1 reverse) + xdg::dblRay *ray; // double precision rays GPRTPrimitiveRef* primitive_refs; int num_faces; // Number of faces in the geometry }; struct dblRayGenData { - dblRay *ray; - dblHit *hit; + xdg::dblRay *ray; + xdg::dblHit *hit; + SurfaceAccelerationStructure* meshid_to_accel_address; // MeshID->TLAS address table to recover volume to trace against }; /* A small structure of constants that can change every frame without rebuilding the shader binding table. (must be 128 bytes or less) */ - -struct dblRayFirePushConstants { +struct dblPushConstants { double tMax; double tMin; + xdg::HitOrientation hitOrientation; }; + +#endif \ No newline at end of file diff --git a/include/xdg/ray_tracing_interface.h b/include/xdg/ray_tracing_interface.h index 9d938978..30a53f17 100644 --- a/include/xdg/ray_tracing_interface.h +++ b/include/xdg/ray_tracing_interface.h @@ -4,16 +4,61 @@ #include #include #include +#include +#include "xdg/error.h" #include "xdg/constants.h" #include "xdg/embree_interface.h" #include "xdg/mesh_manager_interface.h" #include "xdg/primitive_ref.h" #include "xdg/geometry_data.h" + namespace xdg { +struct dblHit; // forward declaration for dblHit + +/** + * @brief Device ray/hit buffer descriptor + * + * This structure provides access to device-allocated ray and hit buffers + * in a backend-agnostic way. The buffers contain XDG's standard ray and hit + * data structures (dblRay and dblHit), regardless of which compute backend + * is being used. + * + * Key design principle: + * - Device pointers are opaque (void*) + * - The data layout is always the XDG types dblRay and dblHit + * - Downstream code can write to these buffers (hopefully) using any compute API + * + * For type-safe access in downstream code: + * - Cast rayDevPtr to (dblRay*) when using C++ or kernels + * - Cast hitDevPtr to (dblHit*) when reading hit results + */ +struct DeviceRayHitBuffers { + void* rayDevPtr; + void* hitDevPtr; + size_t capacity; // Number of rays the buffer can hold + size_t rayStride; // Bytes between ray elements - currently set to sizeof(dblRay) but in theory allows for future flexibility + size_t hitStride; // Bytes between hit elements - currently set to sizeof(dblHit) but in theory allows for future flexibility +}; + +/** + * @brief Callback alias for external ray population + * + * Allows downstream applications to populate ray buffers using their own compute backend + * (GPRT, CUDA, OpenMP) without XDG needing to know the specifics. + * + * The callback receives opaque device pointers and should interpret them according to + * the buffer metadata (stride information). Alternatively, users can rely on the standard + * dblRay/dblHit layouts if they don't need custom padding/alignment. + * + * @param buffer Device ray buffer descriptor with opaque pointers and metadata + * @param numRays Number of rays to generate/populate + */ +using RayPopulationCallback = std::function; + class RayTracer { public: // Constructors/Destructors @@ -73,12 +118,41 @@ class RayTracer { */ virtual void create_global_element_tree() = 0; - // Query Methods + /** + * @brief Check whether a point lies in a specified volume + * + * This method performs a check to see whether a given point is inside a volume provided. + * It computes this by firing a ray from the point and checking whether or not the ray is Entering or Exiting + * the volume boundary. If no direction is provided, a default direction will be used. + * Note - zero length direction vectors are not explicitly checked for internally and should be avoided to avoid causing undefined behavior. + * + * @param[in] tree The TreeID of the volume we are querying against + * @param[in] point The point to be queried + * @param[in] direction (optional) direction to launch a ray in a specified direction - must be non-zero length + * @param[in] exclude_primitives (optional) vector of surface element MeshIDs to exclude from intersection tests + * @return Boolean result of point in volume check + */ virtual bool point_in_volume(TreeID tree, const Position& point, const Direction* direction = nullptr, const std::vector* exclude_primitives = nullptr) const = 0; - + + /** + * @brief Fire a ray against a given volume and return the first hit + * + * This method fires a ray from a given origin in a specified direction against the surfaces of a volume. + * It returns the distance to the closest hit and the MeshID of the surface hit. The user can specify + * a distance limit and whether Entering/Exiting hits should be rejected. + * Note - zero length direction vectors are not explicitly checked for internally and should be avoided to avoid causing undefined behavior. + * + * @param[in] tree The TreeID of the volume we are querying against + * @param[in] origin An array of Position objects representing the starting points of the rays + * @param[in] direction (optional) Direction object to launch a ray in a specified direction + * @param[in] dist_limit (optional) maximum distance to consider for intersections + * @param[in] orientation (optional) flag to consider whether Entering/Exiting hits should be rejected. Defaults to EXITING + * @param[in] exclude_primitives (optional) vector of surface element MeshIDs to exclude from intersection tests + * @return A pair containing the distance to the closest hit and the MeshID of the surface hit + */ virtual std::pair ray_fire(TreeID tree, const Position& origin, const Direction& direction, @@ -123,6 +197,61 @@ class RayTracer { int num_registered_surface_trees() const { return surface_trees_.size(); }; int num_registered_element_trees() const { return element_trees_.size(); }; + + // GPU Ray Tracing Support + + virtual void ray_fire_prepared(const size_t num_rays, + const double dist_limit = INFTY, + HitOrientation orientation = HitOrientation::EXITING) + { + fatal_error("GPU ray tracing not supported with this RayTracer backend"); + } + + virtual void point_in_volume_prepared(const size_t num_points) + { + fatal_error("GPU ray tracing not supported with this RayTracer backend"); + } + + /** + * @brief Check whether the current ray buffer capacity is sufficient for the number of rays requested + * @param[in] num_rays The number of rays to be processed + */ + virtual void check_rayhit_buffer_capacity(const size_t num_rays) { + fatal_error("GPU ray tracing not supported with this RayTracer backend"); + } + + /** + * @brief return device pointers to ray and hit buffers for GPU ray tracing + * @return DeviceRayHitBuffers struct containing device pointers to ray and hit buffers + */ + virtual DeviceRayHitBuffers get_device_rayhit_buffers(const size_t num_rays) { + fatal_error("GPU ray tracing not supported with this RayTracer backend"); + return {}; + } + + /** + * @brief Allocate device ray buffers and populate them via a user-provided callback + * + * This method allows downstream applications to populate ray buffers using any compute + * backend (GPRT, CUDA, HIP, OpenCL, etc.) without coupling them to XDG's internals. + * + * The workflow: + * 1. XDG allocates device memory for rays (if not already large enough) + * 2. XDG passes device pointers to the callback + * 3. User's callback populates the buffers using their preferred compute kernel/shader + * 4. User's callback returns (XDG assumes buffers are now populated) + * + * This avoids unnecessary host-device transfers by allowing users to write directly + * to XDG's device buffers without any host-side transfers. + * + * @param numRays Number of rays to allocate space for + * @param callback Function that will populate the ray buffer. Receives the allocated buffer and ray count. + */ + virtual void populate_rays_external(size_t numRays, + const RayPopulationCallback& callback) { + fatal_error("GPU ray tracing not supported with this RayTracer backend"); + } + protected: // Common functions across RayTracers const double bounding_box_bump(const std::shared_ptr mesh_manager, MeshID volume_id); // return a bump value based on the size of a bounding box (minimum 1e-3). Should this be a part of mesh_manager? @@ -146,8 +275,6 @@ class RayTracer { ElementTreeID next_element_tree_id_ {0}; double numerical_precision_ {1e-3}; }; - } // namespace xdg - #endif // include guard \ No newline at end of file diff --git a/include/xdg/xdg.h b/include/xdg/xdg.h index 82c1f8f2..139f240d 100644 --- a/include/xdg/xdg.h +++ b/include/xdg/xdg.h @@ -10,6 +10,8 @@ namespace xdg { +struct DeviceRayHitBuffers; // forward declaration +struct dblHit; // forward declaration class XDG { public: @@ -63,11 +65,40 @@ next_element(MeshID current_element, const Position& r, const Direction& u) const; +/** + * @brief Check whether a point lies in a specified volume + * + * This method performs a check to see whether a given point is inside a volume provided. + * It computes this by firing a ray from the point and checking whether or not the ray is Entering or Exiting + * the volume boundary. If no direction is provided, a default direction will be used. + * Note - zero length direction vectors are not explicitly checked for internally and should be avoided to avoid causing undefined behavior. + * + * @param[in] tree The TreeID of the volume we are querying against + * @param[in] point The point to be queried + * @param[in] direction (optional) direction to launch a ray in a specified direction - must be non-zero length + * @param[in] exclude_primitives (optional) vector of surface element MeshIDs to exclude from intersection tests + * @return Boolean result of point in volume check + */ bool point_in_volume(MeshID volume, const Position point, const Direction* direction = nullptr, const std::vector* exclude_primitives = nullptr) const; +/** + * @brief Fire a ray against a given volume and return the first hit + * + * This method fires a ray from a given origin in a specified direction against the surfaces of a volume. + * It returns the distance to the closest hit and the MeshID of the surface hit. The user can specify + * a distance limit and whether Entering/Exiting hits should be rejected. + * + * @param[in] volume The MeshID of the volume we are querying against + * @param[in] origin Origin of the ray to be fired + * @param[in] direction (optional) Direction object to launch a ray in a specified direction + * @param[in] dist_limit (optional) maximum distance to consider for intersections + * @param[in] orientation (optional) flag to consider whether Entering/Exiting hits should be rejected. Defaults to EXITING + * @param[in] exclude_primitives (optional) vector of surface element MeshIDs to exclude from intersection tests + * @return A pair containing the distance to the closest hit and the MeshID of the surface hit + */ std::pair ray_fire(MeshID volume, const Position& origin, const Direction& direction, @@ -75,6 +106,33 @@ std::pair ray_fire(MeshID volume, HitOrientation orientation = HitOrientation::EXITING, std::vector* const exclude_primitives = nullptr) const; + +/** + * @brief Call ray fire on pre-populated ray buffers + * + * This method performs a set of ray fire queries on a set of rays that have already been populated on device + * via the external ray population callback method. With GPRT ray tracing this launches the RT pipeline with the number of rays provided. + * + * @param[in] num_rays The number of rays to be processed in the batch + * @param[in] dist_limit (optional) maximum distance to consider for intersections + * @param[in] orientation (optional) flag to consider whether Entering/Exiting hits should be rejected. Defaults to EXITING + * @return Void. Outputs stored in dblHit buffer on device. + */ +void ray_fire_prepared(const size_t num_rays, + const double dist_limit = INFTY, + HitOrientation orientation = HitOrientation::EXITING); + +/** + * @brief Call point_in_volume on pre-populated ray buffers + * + * This method performs a set of point_in_volume queries on a set of points that have already been populated on device + * via the external ray population callback method. With GPRT ray tracing this launches the RT pipeline with the number of points provided. + * + * @param[in] num_points The number of points to be processed in the batch + * @return Void. Outputs stored in dblHit buffer on device. + */ +void point_in_volume_prepared(const size_t num_points); + std::pair closest(MeshID volume, const Position& origin) const; @@ -105,6 +163,18 @@ Direction surface_normal(MeshID surface, ray_tracing_interface_ = ray_tracing_interface; } + // Resize buffers (if necessary) and return device pointers for ray and hit data + DeviceRayHitBuffers get_device_rayhit_buffers(const size_t requiredCapacity) + { + return ray_tracing_interface()->get_device_rayhit_buffers(requiredCapacity); + } + + void populate_rays_external(size_t numRays, + const RayPopulationCallback& callback) + { + return ray_tracing_interface()->populate_rays_external(numRays, callback); + } + // Accessors const std::shared_ptr& ray_tracing_interface() const { return ray_tracing_interface_; @@ -113,6 +183,7 @@ Direction surface_normal(MeshID surface, const std::shared_ptr& mesh_manager() const { return mesh_manager_; } + // Private methods private: double _triangle_volume_contribution(const PrimitiveRef& triangle) const; diff --git a/src/gprt/dbl_deviceCode.slang b/src/gprt/dbl_deviceCode.slang index ff082a45..360d19b5 100644 --- a/src/gprt/dbl_deviceCode.slang +++ b/src/gprt/dbl_deviceCode.slang @@ -1,7 +1,7 @@ #include "../../include/xdg/gprt/shared_structs.h" [[vk::push_constant]] -dblRayFirePushConstants PC; +dblPushConstants PC; struct RayFirePayload { double distance; // Distance to intersection @@ -53,22 +53,27 @@ void ray_fire_miss(inout RayFirePayload payload) { void ray_fire(uniform dblRayGenData record, uniform DPTriangleGeomData mesh) { RayFirePayload payload; uint rayID = DispatchRaysIndex().x; + xdg::dblRay ray = record.ray[rayID]; // recover ray from buffer // Trace the ray into the scene RayDesc rayDesc; - rayDesc.Origin = float3(record.ray[rayID].origin); - rayDesc.Direction = normalize(float3(record.ray[rayID].direction)); - rayDesc.TMin = float(record.ray[rayID].tMin); - rayDesc.TMax = float(record.ray[rayID].tMax); + rayDesc.Origin = float3(ray.origin); + rayDesc.Direction = normalize(float3(ray.direction)); + rayDesc.TMin = float(PC.tMin); + rayDesc.TMax = float(PC.tMax); - SurfaceAccelerationStructure world = record.ray[rayID].volume_accel; + int mesh_id = ray.volume_mesh_id; + // Recover the TLAS we are tracing against for this ray + SurfaceAccelerationStructure world = record.meshid_to_accel_address[mesh_id]; - // Pass the ray's origin and direction to the payload + // Set payload default values payload.distance = -1.0f; payload.surf_id = -1; payload.tlas = world; - TraceRay(world, RAY_FLAG_NONE, 0xff, 0, 1, rayDesc, payload); + if (ray.enabled == 1u) { // skip traversal for rays that are marked as disabled + TraceRay(world, RAY_FLAG_NONE, 0xff, 0, 1, rayDesc, payload); + } // Store the distance to the hit point and the surface ID in buffers for CPU record.hit[rayID].distance = payload.distance; @@ -80,22 +85,26 @@ void ray_fire(uniform dblRayGenData record, uniform DPTriangleGeomData mesh) { void point_in_volume(uniform dblRayGenData record, uniform DPTriangleGeomData mesh) { RayFirePayload payload; uint rayID = DispatchRaysIndex().x; + xdg::dblRay ray = record.ray[rayID]; // recover ray from buffer // Trace the ray into the scene RayDesc rayDesc; - rayDesc.Origin = float3(record.ray[rayID].origin); - rayDesc.Direction = float3(normalize(record.ray[rayID].direction)); - rayDesc.TMin = float(record.ray[rayID].tMin); - rayDesc.TMax = float(record.ray[rayID].tMax); - - SurfaceAccelerationStructure world = record.ray[rayID].volume_accel; - + rayDesc.Origin = float3(ray.origin); + rayDesc.Direction = float3(normalize(ray.direction)); + rayDesc.TMin = float(PC.tMin); + rayDesc.TMax = float(PC.tMax); + + int mesh_id = ray.volume_mesh_id; + // Recover the TLAS we are tracing against for this ray + SurfaceAccelerationStructure world = record.meshid_to_accel_address[mesh_id]; // Pass the ray's origin and direction to the payload payload.surf_id = -1; payload.tlas = world; payload.piv = xdg::PointInVolume::OUTSIDE; // Initialize point in volume check result to outside (0) - TraceRay(world, RAY_FLAG_NONE, 0xff, 0, 1, rayDesc, payload); + if (ray.enabled == 1u) { // skip traversal for rays that are marked as disabled + TraceRay(world, RAY_FLAG_NONE, 0xff, 0, 1, rayDesc, payload); + } record.hit[rayID].surf_id = payload.surf_id; record.hit[rayID].piv = payload.piv; // Point in volume check result @@ -104,10 +113,14 @@ void point_in_volume(uniform dblRayGenData record, uniform DPTriangleGeomData me // ------------------------------------------------- Compute Shaders ------------------------------------------------- /* A shader to compute and store AABB min/maxes in single precision using double precision coords*/ [shader("compute")] -[numthreads(1, 1, 1)] -void -populate_aabbs(uint3 DispatchThreadID: SV_DispatchThreadID, uniform DPTriangleGeomData record) { +[numthreads(64, 1, 1)] +void populate_aabbs(uint3 DispatchThreadID: SV_DispatchThreadID, uniform DPTriangleGeomData record) { int primID = DispatchThreadID.x; + + // guard against more threads than primitives + if (primID >= record.num_faces) + return; + int3 indices = record.index[primID]; double3 A = record.vertex[indices[0]]; double3 B = record.vertex[indices[1]]; @@ -123,7 +136,6 @@ populate_aabbs(uint3 DispatchThreadID: SV_DispatchThreadID, uniform DPTriangleGe // ------------------------------------------------ CUSTOM INTERSECTION SHADERS ------------------------------------------------ - /* 1D ray generation intersection with a double precision triangle using the Plucker intersection algorithm*/ [shader("intersection")] void DPTrianglePluckerIntersection(uniform DPTriangleGeomData record) @@ -134,6 +146,7 @@ void DPTrianglePluckerIntersection(uniform DPTriangleGeomData record) uint rayID = DispatchRaysIndex().x; uint nRays = DispatchRaysDimensions().x; uint flags = RayFlags(); + xdg::dblRay ray = record.ray[rayID]; // recover ray from buffer if (rayID >= nRays) { return; @@ -155,11 +168,11 @@ void DPTrianglePluckerIntersection(uniform DPTriangleGeomData record) double3 v1 = record.vertex[indices[1]]; double3 v2 = record.vertex[indices[2]]; - double3 origin = record.ray[rayID].origin; - double3 direction = record.ray[rayID].direction; + double3 origin = ray.origin; + double3 direction = ray.direction; - double tMin = record.ray[rayID].tMin; - double tMax = record.ray[rayID].tMax; + double tMin = PC.tMin; + double tMax = PC.tMax; const double3 raya = direction; const double3 rayb = cross(direction, origin); @@ -228,8 +241,9 @@ void DPTrianglePluckerIntersection(uniform DPTriangleGeomData record) double3 norm = record.normals[primID]; // recover double precision normal. TODO - Should we calculate from vertices instead? - // sense adjustment of normal - if (record.ray[rayID].volume_tree == record.reverse_tree) + // sense adjustment of normal (per MeshID) + int mesh_sense = record.meshid_to_sense[ray.volume_mesh_id]; // +1 forward, -1 reverse + if (mesh_sense < 0) { norm = -norm; } @@ -238,16 +252,16 @@ void DPTrianglePluckerIntersection(uniform DPTriangleGeomData record) uint hit_kind = norm_dot_dir < 0 ? HIT_KIND_TRIANGLE_FRONT_FACE : HIT_KIND_TRIANGLE_BACK_FACE; - xdg::HitOrientation hitOrientation = record.ray[rayID].hitOrientation; + xdg::HitOrientation hitOrientation = PC.hitOrientation; if (orientation_cull(direction, norm, hitOrientation)) { return; } - for (int i = 0; i < record.ray[rayID].exclude_count; ++i) + for (int i = 0; i < ray.exclude_count; ++i) { - if (record.ray[rayID].exclude_primitives[i] == global_prim_id) { + if (ray.exclude_primitives[i] == global_prim_id) { return; } } @@ -312,4 +326,4 @@ float next_after(float a) { a_++; } return asfloat(a_); -} \ No newline at end of file +} diff --git a/src/gprt/ray_tracer.cpp b/src/gprt/ray_tracer.cpp index bb3ece22..d5567544 100644 --- a/src/gprt/ray_tracer.cpp +++ b/src/gprt/ray_tracer.cpp @@ -1,25 +1,29 @@ #include "xdg/gprt/ray_tracer.h" #include "gprt/gprt.h" - namespace xdg { GPRTRayTracer::GPRTRayTracer() { + gprtRequestRayTypeCount(numRayTypes_); // Set the number of shaders which can be set to the same geometry context_ = gprtContextCreate(); module_ = gprtModuleCreate(context_, dbl_deviceCode); - rayHitBuffers_.capacity = 1; // Preallocate space for 1 ray - rayHitBuffers_.ray = gprtDeviceBufferCreate(context_, rayHitBuffers_.capacity); - rayHitBuffers_.hit = gprtDeviceBufferCreate(context_, rayHitBuffers_.capacity); - rayHitBuffers_.devRayAddr = gprtBufferGetDevicePointer(rayHitBuffers_.ray); - rayHitBuffers_.devHitAddr = gprtBufferGetDevicePointer(rayHitBuffers_.hit); + // Buffer setup + rayHitBuffers_.view.capacity = 1e6; // Preallocate space for 1m rays + rayHitBuffers_.ray = gprtDeviceBufferCreate(context_, rayHitBuffers_.view.capacity); + rayHitBuffers_.hit = gprtDeviceBufferCreate(context_, rayHitBuffers_.view.capacity); + rayHitBuffers_.view.rayDevPtr = gprtBufferGetDevicePointer(rayHitBuffers_.ray); + rayHitBuffers_.view.hitDevPtr = gprtBufferGetDevicePointer(rayHitBuffers_.hit); + rayHitBuffers_.view.rayStride = sizeof(dblRay); + rayHitBuffers_.view.hitStride = sizeof(dblHit); excludePrimitivesBuffer_ = gprtDeviceBufferCreate(context_); // initialise buffer of size 1 + tlas_handle_buffer_ = gprtDeviceBufferCreate(context_); // initialise buffer of size 1 + setup_shaders(); - // Bind the buffers to the RayGenData structure dblRayGenData* rayGenData = gprtRayGenGetParameters(rayGenPrograms_.at(RayGenType::RAY_FIRE)); rayGenData->ray = gprtBufferGetDevicePointer(rayHitBuffers_.ray); @@ -40,17 +44,11 @@ GPRTRayTracer::~GPRTRayTracer() gprtGraphicsSynchronize(context_); gprtComputeSynchronize(context_); - // Destroy TLAS structures for (const auto& [tree, accel] : surface_volume_tree_to_accel_map) { gprtAccelDestroy(accel); } - // Destroy BLAS structures - for (const auto& blas : blas_handles_) { - gprtAccelDestroy(blas); - } - // Destroy Geoms and Types for (const auto& [surf, geom] : surface_to_geometry_map_) { gprtGeomDestroy(geom); @@ -61,6 +59,7 @@ GPRTRayTracer::~GPRTRayTracer() gprtBufferDestroy(rayHitBuffers_.ray); gprtBufferDestroy(rayHitBuffers_.hit); gprtBufferDestroy(excludePrimitivesBuffer_); + gprtBufferDestroy(tlas_handle_buffer_); // Destroy module and context gprtModuleDestroy(module_); @@ -85,9 +84,13 @@ void GPRTRayTracer::setup_shaders() void GPRTRayTracer::init() { + update_tlas_table_(); + // Build the shader binding table (SBT) after all shader programs and acceleration structures are set up gprtBuildShaderBindingTable(context_, GPRT_SBT_ALL); // Note that should we need to update any shaders or acceleration structures, we must rebuild the SBT + + initialized_ = true; } std::pair @@ -164,13 +167,14 @@ GPRTRayTracer::create_surface_tree(const std::shared_ptr& mesh_mana geom_data->normals = gprtBufferGetDevicePointer(normal_buffer); geom_data->primitive_refs = gprtBufferGetDevicePointer(primitive_refs_buffer); geom_data->num_faces = num_faces; - - gprtComputeLaunch(aabbPopulationProgram_, {num_faces, 1, 1}, {1, 1, 1}, *geom_data); + // meshid_to_sense pointer is set after meshid_to_sense_buffer_ is created - GPRTAccel blas = gprtAABBAccelCreate(context_, triangleGeom, buildParams_.buildMode); + constexpr uint32_t threadsPerGroup = 64; // must match [numthreads(64,1,1)] + uint32_t numGroupsX = (num_faces + threadsPerGroup - 1) / threadsPerGroup; + gprtComputeLaunch(aabbPopulationProgram_, {numGroupsX, 1, 1}, {threadsPerGroup, 1, 1}, *geom_data); + GPRTAccel blas = gprtAABBAccelCreate(context_, triangleGeom, buildParams_.buildMode); gprtAccelBuild(context_, blas, buildParams_); - gprt::Instance instance; instance = gprtAccelGetInstance(blas); // create instance of BLAS to be added to TLAS instance.mask = 0xff; // mask can be used to filter instances during ray traversal. 0xff ensures no filtering @@ -184,14 +188,15 @@ GPRTRayTracer::create_surface_tree(const std::shared_ptr& mesh_mana surfaceBlasInstances.push_back(instance); globalBlasInstances_.push_back(instance); - // Always update per-volume info + // Ensure MeshID->sense lookup has an entry for this volume and record its + // orientation sign (+1 forward, -1 reverse) for normal flipping in shader. auto [forward_parent, reverse_parent] = mesh_manager->get_parent_volumes(surf); if (volume_id == forward_parent) { - geom_data->forward_vol = forward_parent; - geom_data->forward_tree = tree; + meshid_to_sense_.resize(static_cast(forward_parent) + 1, 1); + meshid_to_sense_[forward_parent] = 1; } else if (volume_id == reverse_parent) { - geom_data->reverse_vol = reverse_parent; - geom_data->reverse_tree = tree; + meshid_to_sense_.resize(static_cast(reverse_parent) + 1, 1); + meshid_to_sense_[reverse_parent] = -1; } else { fatal_error("Volume {} is not a parent of surface {}", volume_id, surf); } @@ -202,7 +207,19 @@ GPRTRayTracer::create_surface_tree(const std::shared_ptr& mesh_mana GPRTAccel volume_tlas = gprtInstanceAccelCreate(context_, surfaceBlasInstances.size(), instanceBuffer); gprtAccelBuild(context_, volume_tlas, buildParams_); surface_volume_tree_to_accel_map[tree] = volume_tlas; - + surface_tree_to_volume_map_[tree] = volume_id; + if (volume_id >= tlas_handles_.size()) { + tlas_handles_.resize(volume_id + 1, SurfaceAccelerationStructure{}); + } + tlas_handles_[volume_id] = gprtAccelGetDeviceAddress(volume_tlas); + + update_meshid_to_sense_(); + + if (initialized_) { + update_tlas_table_(); + gprtBuildShaderBindingTable(context_, GPRT_SBT_ALL); + } + return tree; } @@ -218,23 +235,23 @@ bool GPRTRayTracer::point_in_volume(SurfaceTreeID tree, const Direction* direction, const std::vector* exclude_primitives) const { - GPRTAccel volume = surface_volume_tree_to_accel_map.at(tree); + MeshID volume = surface_tree_to_volume_map_.at(tree); // recover MeshID of volume to return GPRTAccel on device auto rayGen = rayGenPrograms_.at(RayGenType::POINT_IN_VOLUME); dblRayGenData* rayGenPIVData = gprtRayGenGetParameters(rayGen); + const Direction defaultDir = Direction{1. / std::sqrt(2.0), 1. / std::sqrt(2.0), 0.0}; + // Use provided direction or if Direction == nulptr use default direction Direction directionUsed = (direction != nullptr) ? Direction{direction->x, direction->y, direction->z} - : Direction{1. / std::sqrt(2.0), 1. / std::sqrt(2.0), 0.0}; + : defaultDir; - gprtBufferMap(rayHitBuffers_.ray); // Update the ray input buffer + // Host -> Device buffer mapping/population of raydata for raygen shader + gprtBufferMap(rayHitBuffers_.ray); dblRay* ray = gprtBufferGetHostPointer(rayHitBuffers_.ray); - ray[0].volume_accel = gprtAccelGetDeviceAddress(volume); ray[0].origin = {point.x, point.y, point.z}; ray[0].direction = {directionUsed.x, directionUsed.y, directionUsed.z}; - ray[0].tMax = INFTY; // Set a large distance limit - ray[0].tMin = 0.0; - ray[0].volume_tree = tree; // Set the TreeID of the volume being queried - ray[0].hitOrientation = HitOrientation::ANY; // No orientation culling for point-in-volume check + ray[0].volume_mesh_id = volume; + ray[0].enabled = 1; // Ensure the ray is enabled if (exclude_primitives) { if (!exclude_primitives->empty()) gprtBufferResize(context_, excludePrimitivesBuffer_, exclude_primitives->size(), false); @@ -252,10 +269,15 @@ bool GPRTRayTracer::point_in_volume(SurfaceTreeID tree, } gprtBufferUnmap(rayHitBuffers_.ray); // required to sync buffer back on GPU? - gprtRayGenLaunch1D(context_, rayGen, 1); // Launch raygen shader (entry point to RT pipeline) + dblPushConstants pushConstants; + pushConstants.hitOrientation = HitOrientation::ANY; + pushConstants.tMax = INFTY; + pushConstants.tMin = 0.0; + + gprtRayGenLaunch1D(context_, rayGen, 1, pushConstants); // Launch raygen shader (entry point to RT pipeline) gprtGraphicsSynchronize(context_); // Ensure all GPU operations are complete before returning control flow to CPU - // Retrieve the hit from the dblHit buffer + // Device -> Host buffer mapping to retrieve hit result gprtBufferMap(rayHitBuffers_.hit); dblHit* hit = gprtBufferGetHostPointer(rayHitBuffers_.hit); auto surface = hit[0].surf_id; @@ -279,19 +301,17 @@ std::pair GPRTRayTracer::ray_fire(SurfaceTreeID tree, HitOrientation orientation, std::vector* const exclude_primitives) { - GPRTAccel volume = surface_volume_tree_to_accel_map.at(tree); + MeshID volume = surface_tree_to_volume_map_.at(tree); // recover MeshID of volume to return GPRTAccel on device auto rayGen = rayGenPrograms_.at(RayGenType::RAY_FIRE); dblRayGenData* rayGenData = gprtRayGenGetParameters(rayGen); - - gprtBufferMap(rayHitBuffers_.ray); // Update the ray input buffer + + // Host -> Device buffer mapping/population of raydata for raygen shader + gprtBufferMap(rayHitBuffers_.ray); dblRay* ray = gprtBufferGetHostPointer(rayHitBuffers_.ray); - ray[0].volume_accel = gprtAccelGetDeviceAddress(volume); ray[0].origin = {origin.x, origin.y, origin.z}; ray[0].direction = {direction.x, direction.y, direction.z}; - ray[0].tMax = dist_limit; - ray[0].tMin = 0.0; - ray[0].hitOrientation = orientation; // Set orientation for the ray - ray[0].volume_tree = tree; // Set the TreeID of the volume being queried + ray[0].volume_mesh_id = volume; + ray[0].enabled = 1; // Ensure the ray is enabled if (exclude_primitives) { if (!exclude_primitives->empty()) gprtBufferResize(context_, excludePrimitivesBuffer_, exclude_primitives->size(), false); @@ -309,10 +329,16 @@ std::pair GPRTRayTracer::ray_fire(SurfaceTreeID tree, } gprtBufferUnmap(rayHitBuffers_.ray); // required to sync buffer back on GPU? - gprtRayGenLaunch1D(context_, rayGen, 1); // Launch raygen shader (entry point to RT pipeline) + // Set push constants (same for every ray) + dblPushConstants pushConstants; + pushConstants.hitOrientation = orientation; + pushConstants.tMax = dist_limit; + pushConstants.tMin = 0.0; + + gprtRayGenLaunch1D(context_, rayGen, 1, pushConstants); // Launch raygen shader (entry point to RT pipeline) gprtGraphicsSynchronize(context_); // Ensure all GPU operations are complete before returning control flow to CPU - // Retrieve the hit from the dblHit buffer + // Device -> Host buffer mapping to retrieve hit result gprtBufferMap(rayHitBuffers_.hit); dblHit* hit = gprtBufferGetHostPointer(rayHitBuffers_.hit); auto distance = hit[0].distance; @@ -326,7 +352,46 @@ std::pair GPRTRayTracer::ray_fire(SurfaceTreeID tree, if (exclude_primitives) exclude_primitives->push_back(primitive_id); return {distance, surface}; } - + +void +GPRTRayTracer::ray_fire_prepared(const size_t num_rays, + const double dist_limit, + HitOrientation orientation) +{ + if (num_rays == 0) return; // no work to do. Early exit + + check_rayhit_buffer_capacity(num_rays); + + auto rayGen = rayGenPrograms_.at(RayGenType::RAY_FIRE); + + dblPushConstants pushConstants; + pushConstants.tMax = dist_limit; + pushConstants.tMin = 0.0; + pushConstants.hitOrientation = orientation; // Set orientation for the ray + + gprtRayGenLaunch1D(context_, rayGen, num_rays, pushConstants); + gprtGraphicsSynchronize(context_); + return; +} + +void +GPRTRayTracer::point_in_volume_prepared(const size_t num_points) +{ + if (num_points == 0) return; // no work to do. Early exit + + check_rayhit_buffer_capacity(num_points); + auto rayGen = rayGenPrograms_.at(RayGenType::POINT_IN_VOLUME); + + dblPushConstants pushConstants; + pushConstants.tMax = INFTY; + pushConstants.tMin = 0.0; + pushConstants.hitOrientation = HitOrientation::ANY; // Set orientation for the ray + + gprtRayGenLaunch1D(context_, rayGen, num_points, pushConstants); + gprtGraphicsSynchronize(context_); + return; +} + void GPRTRayTracer::create_global_surface_tree() { // Create a TLAS (Top-Level Acceleration Structure) for all the volumes @@ -341,16 +406,22 @@ void GPRTRayTracer::create_global_surface_tree() global_surface_accel_ = global_accel; } -void GPRTRayTracer::check_ray_buffer_capacity(size_t N) +void GPRTRayTracer::check_rayhit_buffer_capacity(const size_t N) { - if (N <= rayHitBuffers_.capacity) return; // current capacity is sufficient + if (N <= rayHitBuffers_.view.capacity) return; // current capacity is sufficient - // Resize buffers to accommodate N rays - size_t newCapacity = std::max(N, rayHitBuffers_.capacity * 2); // double the capacity or set to N, whichever is larger + // Resize buffers to accommodate N rays - double the capacity or set to N, whichever is larger + size_t newCapacity = std::max(N, rayHitBuffers_.view.capacity * 2); gprtBufferResize(context_, rayHitBuffers_.ray, newCapacity, false); gprtBufferResize(context_, rayHitBuffers_.hit, newCapacity, false); - rayHitBuffers_.capacity = newCapacity; + rayHitBuffers_.view.capacity = newCapacity; + + // Get fresh device pointers after resize + rayHitBuffers_.view.rayDevPtr = gprtBufferGetDevicePointer(rayHitBuffers_.ray); + rayHitBuffers_.view.hitDevPtr = gprtBufferGetDevicePointer(rayHitBuffers_.hit); + rayHitBuffers_.view.rayStride = sizeof(dblRay); + rayHitBuffers_.view.hitStride = sizeof(dblHit); // Since we have resized the ray buffers, we need to update the geom_data->rayIn pointers in all geometries too for (auto const& [surf, geom] : surface_to_geometry_map_) { @@ -363,9 +434,75 @@ void GPRTRayTracer::check_ray_buffer_capacity(size_t N) dblRayGenData* rayGenData = gprtRayGenGetParameters(rayGen); rayGenData->ray = gprtBufferGetDevicePointer(rayHitBuffers_.ray); rayGenData->hit = gprtBufferGetDevicePointer(rayHitBuffers_.hit); + rayGenData->meshid_to_accel_address = gprtBufferGetDevicePointer(tlas_handle_buffer_); } gprtBuildShaderBindingTable(context_, static_cast(GPRT_SBT_GEOM | GPRT_SBT_RAYGEN)); } +// Update the TLAS table (MeshID -> SurfaceAccelerationStructure) buffer on the device +void GPRTRayTracer::update_tlas_table_() +{ + upload_device_buffer_(tlas_handle_buffer_, tlas_handles_); + + for (auto type : {RayGenType::RAY_FIRE, RayGenType::POINT_IN_VOLUME}) { + auto* raygendata = gprtRayGenGetParameters(rayGenPrograms_.at(type)); + raygendata->meshid_to_accel_address = gprtBufferGetDevicePointer(tlas_handle_buffer_); + } +} + +void GPRTRayTracer::update_meshid_to_sense_() +{ + upload_device_buffer_(meshid_to_sense_buffer_, meshid_to_sense_); + + for (auto const& [surf, geom] : surface_to_geometry_map_) { + DPTriangleGeomData* geom_data = gprtGeomGetParameters(geom); + geom_data->meshid_to_sense = gprtBufferGetDevicePointer(meshid_to_sense_buffer_); + } +} + + +DeviceRayHitBuffers GPRTRayTracer::get_device_rayhit_buffers(const size_t N) +{ + check_rayhit_buffer_capacity(N); + return rayHitBuffers_.view; +} + +void GPRTRayTracer::populate_rays_external(size_t numRays, + const RayPopulationCallback& callback) +{ + if (numRays == 0) { + warning("Warning number of rays passed to populate_rays_external is 0. No work to be done."); + return; + } + + // Ensure device buffers are large enough + check_rayhit_buffer_capacity(numRays); + + // Use the user callback to populate the rays directly on the device + callback(rayHitBuffers_.view, numRays); + + // After callback returns, we assume the ray buffer is populated and ready to trace + // Note: The callback is responsible for synchronization if using an async API +} + +void GPRTRayTracer::transfer_hits_buffer_to_host(const size_t num_rays, + std::vector& hits) +{ + hits.clear(); // Ensure hits vector is empty before populating + if (num_rays == 0) { + warning("Warning number of rays passed to transfer_hits_buffer_to_host is 0. No work to be done."); + return; + } + if (num_rays > rayHitBuffers_.view.capacity) { + fatal_error("Requested {} hits, but hit buffer capacity is {}", num_rays, rayHitBuffers_.view.capacity); + } + + hits.resize(num_rays); + gprtBufferMap(rayHitBuffers_.hit); + dblHit* hit = gprtBufferGetHostPointer(rayHitBuffers_.hit); + std::copy(hit, hit + num_rays, hits.begin()); + gprtBufferUnmap(rayHitBuffers_.hit); +} + } // namespace xdg diff --git a/src/tetrahedron_contain.cpp b/src/tetrahedron_contain.cpp index 743ea1d0..c6f43a2a 100644 --- a/src/tetrahedron_contain.cpp +++ b/src/tetrahedron_contain.cpp @@ -13,7 +13,10 @@ bool plucker_tet_containment_test(const Position& point, const Position& v1, const Position& v2, const Position& v3) { - using namespace linalg::aliases; + // explicit namespace usage to avoid clash with GPRT math types + using linalg::aliases::double3x3; + using linalg::aliases::double3; + using linalg::aliases::double4; // Create matrix T = [v1 - v0, v2 - v0, v3 - v0] Vec3da e0 = v1 - v0; Vec3da e1 = v2 - v0; diff --git a/src/xdg.cpp b/src/xdg.cpp index 371a9a1c..46cb7795 100644 --- a/src/xdg.cpp +++ b/src/xdg.cpp @@ -235,8 +235,22 @@ XDG::ray_fire(MeshID volume, HitOrientation orientation, std::vector* const exclude_primitives) const { - TreeID scene = volume_to_surface_tree_map_.at(volume); - return ray_tracing_interface()->ray_fire(scene, origin, direction, dist_limit, orientation, exclude_primitives); + TreeID tree = volume_to_surface_tree_map_.at(volume); + return ray_tracing_interface()->ray_fire(tree, origin, direction, dist_limit, orientation, exclude_primitives); +} + +void +XDG::ray_fire_prepared(const size_t num_rays, + const double dist_limit, + HitOrientation orientation) +{ + return ray_tracing_interface()->ray_fire_prepared(num_rays, dist_limit, orientation); +} + +void +XDG::point_in_volume_prepared(const size_t num_points) +{ + return ray_tracing_interface()->point_in_volume_prepared(num_points); } std::pair XDG::closest(MeshID volume, diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index d7286555..8166eee9 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -37,6 +37,11 @@ if (XDG_ENABLE_MOAB) list(APPEND TEST_NAMES test_overlap_check) endif() +# This test really should be appended when any GPU library is enabled +if (XDG_ENABLE_GPRT) + list(APPEND TEST_NAMES test_direct_ray_buffer_access) +endif() + foreach(test ${TEST_NAMES}) add_executable(${test} ${test}.cpp) target_link_libraries(${test} xdg Catch2::Catch2WithMain) @@ -48,6 +53,18 @@ foreach(test ${TEST_NAMES}) TEST_PREFIX "${test}::") endforeach() +if (XDG_ENABLE_GPRT) + embed_devicecode( + OUTPUT_TARGET + test_direct_ray_buffer_access_deviceCode + HEADERS + ${CMAKE_CURRENT_SOURCE_DIR}/test_direct_ray_buffer_access_shared.h + SOURCES + ${CMAKE_CURRENT_SOURCE_DIR}/test_direct_ray_buffer_access_deviceCode.slang + ) + target_link_libraries(test_direct_ray_buffer_access test_direct_ray_buffer_access_deviceCode) +endif() + set( TEST_FILES diff --git a/tests/test_direct_ray_buffer_access.cpp b/tests/test_direct_ray_buffer_access.cpp new file mode 100644 index 00000000..8b75900c --- /dev/null +++ b/tests/test_direct_ray_buffer_access.cpp @@ -0,0 +1,201 @@ +// for testing +#include +#include +#include +#include + +// xdg includes +#include "xdg/constants.h" +#include "xdg/mesh_manager_interface.h" +#include "xdg/gprt/ray_tracer.h" +#include "xdg/xdg.h" +#include "mesh_mock.h" +#include "test_direct_ray_buffer_access_shared.h" +#include "util.h" +#include "gprt.h" + +#include + +using namespace xdg; +using namespace xdg::test; + +extern GPRTProgram test_direct_ray_buffer_access_deviceCode; + +static RayPopulationCallback make_populate_callback(const std::vector& origins, + const std::vector& directions, + const std::vector& volume_ids, + GPRTContext context, + GPRTComputeOf packRays) { + return [&origins, &directions, volume_ids, context, packRays] + (const DeviceRayHitBuffers& buffer, size_t num_rays) { + + // When passing arrays to the callback, ensure they are the correct size + assert(origins.size() == num_rays); + assert(directions.size() == num_rays); + + // Convert to double3 for use on GPU + std::vector origins_device(num_rays); + std::vector directions_device(num_rays); + for (size_t i = 0; i < num_rays; ++i) { + origins_device[i] = {origins[i].x, origins[i].y, origins[i].z}; + directions_device[i] = {directions[i].x, directions[i].y, directions[i].z}; + } + + auto origins_buffer = gprtDeviceBufferCreate(context, num_rays, origins_device.data()); + auto directions_buffer = gprtDeviceBufferCreate(context, num_rays, directions_device.data()); + auto volume_ids_buffer = gprtDeviceBufferCreate(context, num_rays, volume_ids.data()); + + constexpr uint32_t threads_per_group = 256; + const uint32_t groups = static_cast((num_rays + threads_per_group - 1) / threads_per_group); + + ExternalRayParams params = {}; + params.xdgRays = static_cast(buffer.rayDevPtr); + params.origins = gprtBufferGetDevicePointer(origins_buffer); + params.directions = gprtBufferGetDevicePointer(directions_buffer); + params.num_rays = static_cast(num_rays); + params.total_threads = groups * threads_per_group; + params.volume_mesh_ids = gprtBufferGetDevicePointer(volume_ids_buffer); // Pass array of volume IDs to compute shader + params.enabled = 1u; + + gprtComputeLaunch(packRays, + { groups, 1, 1 }, + { threads_per_group, 1, 1 }, + params); + gprtComputeSynchronize(context); + + gprtBufferDestroy(origins_buffer); + gprtBufferDestroy(directions_buffer); + if (volume_ids_buffer) { + gprtBufferDestroy(volume_ids_buffer); + } + }; +} + +// This is a GPU only test - skip if no GPU ray tracing backends are enabled +TEMPLATE_TEST_CASE("Ray Fire with external populated rays", "[rayfire][mock]", + GPRT_Raytracer) +{ + constexpr auto rt_backend = TestType::value; + + DYNAMIC_SECTION(fmt::format("Backend = {}", rt_backend)) { + check_ray_tracer_supported(rt_backend); // skip if backend not enabled at configuration time + std::shared_ptr xdg = XDG::create(MeshLibrary::MOAB, rt_backend); + REQUIRE(xdg->ray_tracing_interface()->library() == rt_backend); + REQUIRE(xdg->mesh_manager()->mesh_library() == MeshLibrary::MOAB); + const auto& mesh_manager = xdg->mesh_manager(); + mesh_manager->load_file("jezebel.h5m"); + mesh_manager->init(); + xdg->prepare_raytracer(); + + std::vector origins; + std::vector directions; + size_t N = 64; + make_rays(N, origins, directions); + + auto gprt_rt = std::dynamic_pointer_cast(xdg->ray_tracing_interface()); + REQUIRE(gprt_rt); + + std::vector volumes = mesh_manager->volumes(); + REQUIRE(volumes.size() >= 2); + const MeshID volume_id = volumes[0]; + const MeshID volume_id_alt = volumes[1]; + GPRTContext context = gprt_rt->context(); + GPRTModule module = gprtModuleCreate(context, test_direct_ray_buffer_access_deviceCode); + auto packRays = gprtComputeCreate(context, module, "pack_external_rays"); + + std::vector expected_distances(N, INFTY); + std::vector expected_surfaces(N, ID_NONE); + std::vector volume_ids(N, volume_id); + for (size_t i = 0; i < N; ++i) { + volume_ids[i] = (i % 2 == 0) ? volume_id : volume_id_alt; // Volume IDs alternating between two volumes + const auto [dist, surf] = xdg->ray_fire(volume_ids[i], origins[i], directions[i]); + expected_distances[i] = dist; + expected_surfaces[i] = surf; + } + + // Create callback to populate rays on device + RayPopulationCallback populate_callback = make_populate_callback(origins, + directions, + volume_ids, + context, + packRays); + + // Populate rays via external API + xdg->populate_rays_external(N, populate_callback); + + xdg->ray_fire_prepared(N); + std::vector hits; + gprt_rt->transfer_hits_buffer_to_host(N, hits); + + REQUIRE(hits.size() == N); + for (size_t i = 0; i < N; ++i) { + REQUIRE(hits[i].surf_id == expected_surfaces[i]); + if (expected_surfaces[i] != ID_NONE) { + REQUIRE_THAT(hits[i].distance, Catch::Matchers::WithinAbs(expected_distances[i], 1e-6)); + } + } + + gprtComputeDestroy(packRays); + gprtModuleDestroy(module); + } +} + +TEMPLATE_TEST_CASE("Point-in-volume with external populated rays", "[piv][mock]", + GPRT_Raytracer) +{ + constexpr auto rt_backend = TestType::value; + + DYNAMIC_SECTION(fmt::format("Backend = {}", rt_backend)) { + check_ray_tracer_supported(rt_backend); // skip if backend not enabled at configuration time + std::shared_ptr xdg = XDG::create(MeshLibrary::MOAB, rt_backend); + REQUIRE(xdg->ray_tracing_interface()->library() == rt_backend); + REQUIRE(xdg->mesh_manager()->mesh_library() == MeshLibrary::MOAB); + const auto& mesh_manager = xdg->mesh_manager(); + mesh_manager->load_file("jezebel.h5m"); + mesh_manager->init(); + xdg->prepare_raytracer(); + + std::vector points; + std::vector directions; + size_t N = 64; + make_points(N, points, directions); + + auto gprt_rt = std::dynamic_pointer_cast(xdg->ray_tracing_interface()); + REQUIRE(gprt_rt); + + std::vector volumes = mesh_manager->volumes(); + REQUIRE(volumes.size() >= 2); + const MeshID volume_id = volumes[0]; + const MeshID volume_id_alt = volumes[1]; + GPRTContext context = gprt_rt->context(); + GPRTModule module = gprtModuleCreate(context, test_direct_ray_buffer_access_deviceCode); + auto packRays = gprtComputeCreate(context, module, "pack_external_rays"); + + std::vector expected_piv(N, 0); + std::vector volume_ids(N, volume_id); + for (size_t i = 0; i < N; ++i) { + volume_ids[i] = (i % 2 == 0) ? volume_id : volume_id_alt; + expected_piv[i] = static_cast(xdg->point_in_volume(volume_ids[i], points[i], &directions[i])); + } + + RayPopulationCallback populate_callback = make_populate_callback(points, + directions, + volume_ids, + context, + packRays); + xdg->populate_rays_external(N, populate_callback); + + xdg->point_in_volume_prepared(N); + std::vector hits; + gprt_rt->transfer_hits_buffer_to_host(N, hits); + + REQUIRE(hits.size() == N); + for (size_t i = 0; i < N; ++i) { + const auto expected = expected_piv[i] ? xdg::PointInVolume::INSIDE : xdg::PointInVolume::OUTSIDE; // convert back to enum + REQUIRE(hits[i].piv == expected); + } + + gprtComputeDestroy(packRays); + gprtModuleDestroy(module); + } +} diff --git a/tests/test_direct_ray_buffer_access_deviceCode.slang b/tests/test_direct_ray_buffer_access_deviceCode.slang new file mode 100644 index 00000000..b61e6126 --- /dev/null +++ b/tests/test_direct_ray_buffer_access_deviceCode.slang @@ -0,0 +1,24 @@ +#include "test_direct_ray_buffer_access_shared.h" + +[shader("compute")] +[numthreads(256, 1, 1)] +void pack_external_rays(uint3 DispatchThreadID: SV_DispatchThreadID, + uniform ExternalRayParams extParams) +{ + uint globalThreadID = DispatchThreadID.x; + uint stride = extParams.total_threads; // Groups * 256 + + // Grid-stride loop: each thread handles ray idx, idx+stride, idx+2*stride, ... + for (uint idx = globalThreadID; idx < extParams.num_rays; idx += stride) + { + xdg::dblRay r; + r.origin = extParams.origins[idx]; + r.direction = extParams.directions[idx]; + r.exclude_primitives = nullptr; + r.exclude_count = 0; + r.volume_mesh_id = extParams.volume_mesh_ids[idx]; // Set volume mesh ID per ray + r.enabled = extParams.enabled; + + extParams.xdgRays[idx] = r; // Write to device ray buffer + } +} diff --git a/tests/test_direct_ray_buffer_access_shared.h b/tests/test_direct_ray_buffer_access_shared.h new file mode 100644 index 00000000..e8d5ceb3 --- /dev/null +++ b/tests/test_direct_ray_buffer_access_shared.h @@ -0,0 +1,13 @@ +#include "gprt.h" + +#include "../include/xdg/gprt/ray.h" + +struct ExternalRayParams { + xdg::dblRay* xdgRays; + double3* origins; + double3* directions; + uint num_rays; + uint total_threads; + int32_t* volume_mesh_ids; + uint enabled; +}; diff --git a/tests/test_files b/tests/test_files index a3caf0af..eb0334e7 160000 --- a/tests/test_files +++ b/tests/test_files @@ -1 +1 @@ -Subproject commit a3caf0af3f128944c4d6eac93b481df6e4efd97c +Subproject commit eb0334e7bde416845bc28aeb91b192f44ee35a78 diff --git a/tests/test_point_in_volume.cpp b/tests/test_point_in_volume.cpp index ae34e823..63acb882 100644 --- a/tests/test_point_in_volume.cpp +++ b/tests/test_point_in_volume.cpp @@ -24,7 +24,6 @@ TEMPLATE_TEST_CASE("Point-in-volume on MeshMock", "[piv][mock]", check_ray_tracer_supported(rt_backend); // skip if backend not enabled at configuration time auto rti = create_raytracer(rt_backend); REQUIRE(rti); - rti->init(); // Keep MeshMock usage consistent across backends auto mm = std::make_shared(false); @@ -77,4 +76,4 @@ TEMPLATE_TEST_CASE("Point-in-volume on MeshMock", "[piv][mock]", result = rti->point_in_volume(volume_tree, point, &dir); REQUIRE(result == false); } -} +} \ No newline at end of file diff --git a/tests/test_ray_fire.cpp b/tests/test_ray_fire.cpp index be816c36..a730c4c1 100644 --- a/tests/test_ray_fire.cpp +++ b/tests/test_ray_fire.cpp @@ -4,7 +4,6 @@ #include #include - // xdg includes #include "xdg/constants.h" #include "xdg/mesh_manager_interface.h" diff --git a/tests/util.h b/tests/util.h index 579ff652..58d81d8f 100644 --- a/tests/util.h +++ b/tests/util.h @@ -96,3 +96,32 @@ create_raytracer(xdg::RTLibrary rt) { return nullptr; } + +inline void make_rays(size_t N, std::vector& origins, std::vector& directions) +{ + origins.clear(); + directions.clear(); + origins.reserve(N); + directions.reserve(N); + for (size_t i = 0; i < N; ++i) { + int axis = static_cast(i % 3); + double s = (i % 2) ? 1.0 : -1.0; + origins.push_back({0.0, 0.0, 0.0}); + if (axis == 0) directions.push_back({s, 0.0, 0.0}); + else if (axis == 1) directions.push_back({0.0, s, 0.0}); + else directions.push_back({0.0, 0.0, s}); + } +} + +inline void make_points(size_t N, std::vector& points, std::vector& directions) +{ + points.resize(N); + directions.resize(N); + for (size_t i = 0; i < N; ++i) { + // evens inside (origin), odds just outside +X + points[i] = (i % 2 == 0) ? xdg::Position{0.0, 0.0, 0.0} : xdg::Position{5.1, 0.0, 0.0}; + // alternate ±X directions + directions[i] = (i % 2 == 0) ? xdg::Direction{1.0, 0.0, 0.0} + : xdg::Direction{-1.0, 0.0, 0.0}; + } +} diff --git a/vendor/GPRT b/vendor/GPRT index f1e95e41..405d9ee9 160000 --- a/vendor/GPRT +++ b/vendor/GPRT @@ -1 +1 @@ -Subproject commit f1e95e4188cde591547d6b4a33a70bf2afaeec59 +Subproject commit 405d9ee9f5ee8e1a0455f776f9e2c3adffb64160