pushing to work

util/grid/svo.hpp

@@ -15,6 +15,7 @@
 
 class VoxelData;
 
+constexpr float EPSILON = 0.0000000000000000000000001;
 static const size_t CompressionBlockSize = 64*1024*1024;
 
 class VoxelOctree {
@@ -24,14 +25,74 @@ private:
     std::vector<uint32_t> _octree;
     VoxelData* _voxels;
     Vec3f _center;
+
     size_t buildOctree(ChunkedAllocator<uint32_t>& allocator, int x, int y, int z, int size, size_t descriptorIndex) {
         _voxels->prepateDataAccess(x, y, z, size);
 
         int halfSize = size >> 1;
-        const std::array<int, 8> posX = { x + halfSize, x, x+halfSize, x, x+halfSize, x, x+halfSize, x};
+        const std::array<Vec3f, 8> childPositions = {
-        const std::array<int, 8> posY = { y + halfSize, y+halfSize, y, y, y+halfSize, y+halfSize, y, y};
+            Vec3f{x + halfSize, y + halfSize, z + halfSize},
-        const std::array<int, 8> posZ = { z + halfSize, z+halfSize, z+halfSize, z+halfSize, z, z, z, z};
+            Vec3f{x,           y + halfSize, z + halfSize},
+            Vec3f{x + halfSize, y,           z + halfSize},
+            Vec3f{x,           y,           z + halfSize},
+            Vec3f{x + halfSize, y + halfSize, z},
+            Vec3f{x,           y + halfSize, z},
+            Vec3f{x + halfSize, y,           z},
+            Vec3f{x,           y,           z}
+        };
+        uint64_t childOffset = static_cast<uint64_t>(allocator.size()) - descriptorIndex;
 
+        int childCount = 0;
+        std::array<int, 8> childIndices{};
+        uint32_t childMask = 0;
+
+        for (int i = 0; i < 8; ++i) {
+            if (_voxels->cubeContainsVoxelsDestructive(childPositions[i].x, childPositions[i].y, childPositions[i].z, halfSize)) {
+                childMask |= 128 >> i;
+                childIndices[childCount++] = i;
+            }
+        }
+
+        bool hasLargeChildren = false;
+        uint32_t leafMask;
+        if (halfSize == 1) {
+            leafMask = 0;
+            for (int i = 0; i < childCount; ++i) {
+                int idx = childIndices[childCount - i - 1];
+                allocator.pushBack(_voxels->getVoxelDestructive(childPositions[idx].x, childPositions[idx].y, childPositions[idx].z));
+            }
+        } else {
+            leafMask = childMask;
+            for (int i = 0; i < childCount; ++i) allocator.pushBack(0);
+            std::array<uint64_t, 8> granChildOffsets{};
+            uint64_t delta = 0;
+            uint64_t insertionCount = allocator.insertionCount();
+
+            for (int i = 0; i < childCount; ++i) {
+                int idx = childIndices[childCount - i - 1];
+                granChildOffsets[i] = delta + buildOctree(allocator, childPositions[idx].x, childPositions[idx].y, childPositions[idx].z, halfSize, descriptorIndex + childOffset + i);
+                delta += allocator.insertionCount() - insertionCount;
+                insertionCount = allocator.insertionCount();
+                if (granChildOffsets[i] > 0x3FFF) hasLargeChildren = true;
+            }
+
+            for (int i = 0; i < childCount; ++i) {
+                uint64_t childIdx = descriptorIndex + childOffset + i;
+                uint64_t offset = granChildOffsets[i];
+
+                if (hasLargeChildren) {
+                    offset += childCount - i;
+                    allocator.insert(childIdx + 1, static_cast<uint32_t>(offset));
+                    allocator[childIdx] |= 0x20000;
+                    offset >>= 32;
+                }
+                allocator[childIdx] |= static_cast<uint32_t>(offset << 18);
+            }
+        }
+
+        allocator[descriptorIndex] = (childMask << 8) | leafMask;
+        if (hasLargeChildren) allocator[descriptorIndex] |= 0x10000;
+        return childOffset;
     }
 public:
     VoxelOctree(const std::string& path) : _voxels(nullptr) {
@@ -114,7 +175,35 @@ public:
         }
     }
     
-    bool rayMarch(const Vec3f& origin, const Vec3f& dest, float rayScale, uint32_t& normal, float& t);
+    bool rayMarch(const Vec3f& origin, const Vec3f& dest, float rayScale, uint32_t& normal, float& t) {
+        struct StackEntry {
+            uint64_t offset;
+            float maxT;
+        };
+
+        std::array<StackEntry, MaxScale + 1> rayStack;
+
+        Vec3 invAbsD = -dest.abs().safeInverse();
+        uint8_t octantMask = dest.calculateOctantMask();
+        Vec3f bT = invAbsD * origin;
+        if (dest.x > 0) { bT.x = 3.0f * invAbsD.x - bT.x;}
+        if (dest.y > 0) { bT.y = 3.0f * invAbsD.y - bT.y;}
+        if (dest.z > 0) { bT.z = 3.0f * invAbsD.z - bT.z;}
+
+        float minT = (2.0f * invAbsD - bT).maxComp();
+        float maxT = (invAbsD - bT).minComp();
+        minT = std::max(minT, 0.0f);
+        uint32_t curr = 0;
+        uint64_t par = 0;
+
+        Vec3 pos(1.0f);
+        
+        int idx = 0;
+        Vec3 centerT = 1.5f * invAbsD - bT;
+        if (centerT.x > minT) { idx ^= 1; pos.x = 1.5f; }
+        if (centerT.y > minT) { idx ^= 2; pos.y = 1.5f; }
+        if (centerT.z > minT) { idx ^= 4; pos.z = 1.5f; }
+    }
 
     Vec3f center() const {
         return _center;

util/ray3.hpp

@@ -71,4 +71,6 @@ public:
     }
 };
 
+using Ray3f = Ray3<float>;
+
 #endif

util/vectorlogic/vec3.hpp

@@ -242,20 +242,33 @@ public:
                std::abs(z - other.z) < epsilon;
     }
     
-    friend Vec3 operator+(T scalar, const Vec3& vec) {
+    // Template friend operators to allow different scalar types
-        return Vec3(scalar + vec.x, scalar + vec.y, scalar + vec.z);
+    template<typename S>
+    friend Vec3<T> operator+(S scalar, const Vec3<T>& vec) {
+        return Vec3<T>(static_cast<T>(scalar) + vec.x, 
+                      static_cast<T>(scalar) + vec.y, 
+                      static_cast<T>(scalar) + vec.z);
     }
     
-    friend Vec3 operator-(T scalar, const Vec3& vec) {
+    template<typename S>
-        return Vec3(scalar - vec.x, scalar - vec.y, scalar - vec.z);
+    friend Vec3<T> operator-(S scalar, const Vec3<T>& vec) {
+        return Vec3<T>(static_cast<T>(scalar) - vec.x, 
+                      static_cast<T>(scalar) - vec.y, 
+                      static_cast<T>(scalar) - vec.z);
     }
     
-    friend Vec3 operator*(T scalar, const Vec3& vec) {
+    template<typename S>
-        return Vec3(scalar * vec.x, scalar * vec.y, scalar * vec.z);
+    friend Vec3<T> operator*(S scalar, const Vec3<T>& vec) {
+        return Vec3<T>(static_cast<T>(scalar) * vec.x, 
+                      static_cast<T>(scalar) * vec.y, 
+                      static_cast<T>(scalar) * vec.z);
     }
     
-    friend Vec3 operator/(T scalar, const Vec3& vec) {
+    template<typename S>
-        return Vec3(scalar / vec.x, scalar / vec.y, scalar / vec.z);
+    friend Vec3<T> operator/(S scalar, const Vec3<T>& vec) {
+        return Vec3<T>(static_cast<T>(scalar) / vec.x, 
+                      static_cast<T>(scalar) / vec.y, 
+                      static_cast<T>(scalar) / vec.z);
     }
 
     Vec3 reflect(const Vec3& normal) const {
@@ -336,6 +349,30 @@ public:
         return (&x)[index];
     }
     
+    Vec3 safeInverse(float epsilon = 1e-10f) const {
+        return Vec3(
+            1 / (std::abs(x) < epsilon ? std::copysign(epsilon, x) : x),
+            1 / (std::abs(y) < epsilon ? std::copysign(epsilon, y) : y),
+            1 / (std::abs(z) < epsilon ? std::copysign(epsilon, z) : z)
+        );
+    }
+
+    uint8_t calculateOctantMask() const {
+        uint8_t mask = 0;
+        if (x > 0.0f) mask |= 1;
+        if (y > 0.0f) mask |= 2;
+        if (z > 0.0f) mask |= 4;
+        return mask;
+    }
+    
+    float maxComp() const { 
+        return std::max({x, y, z}); 
+    }
+    
+    float minComp() const { 
+        return std::min({x, y, z}); 
+    }
+
     std::string toString() const {
         return "(" + std::to_string(x) + ", " + std::to_string(y) + ", " + std::to_string(z) + ")";
     }