some memory changes and an attempt at fillplanet

tests/planet.cpp

@@ -41,6 +41,7 @@ private:
     std::chrono::steady_clock::time_point lastStatsUpdate;
     std::string cachedStats;
     bool statsNeedUpdate = true;
+    float framerate = 60.0;
 
     enum class DebugColorMode {
         BASE,
@@ -264,6 +265,7 @@ public:
             ImGui::DragFloat("Movement Speed", &cam.movementSpeed, 0.1f, 1.0f, 500.0f);
             ImGui::DragFloat("Rotation Speed", &cam.rotationSpeed, M_1_PI, M_1_PI, M_PI);
             ImGui::InputFloat("Rotation Distance", &rotationRadius, 10, 100);
+            ImGui::InputFloat("Max Framerate", &framerate, 1, 10);
 
             ImGui::Checkbox("Use Slower Render", &slowRender);
 
@@ -508,11 +510,12 @@ public:
         sim.grid.setLODMinDistance(lodDist);
         sim.grid.setLODFalloff(lodDropoff);
         sim.grid.setMaxDistance(maxViewDistance);
+        float invFrameRate = 1 / framerate;
         
         if (slowRender) {
-            currentPreviewFrame = sim.grid.renderFrame(cam, outHeight, outWidth, frame::colormap::RGB, rayCount, reflectCount, globalIllumination, useLod);
+            currentPreviewFrame = sim.grid.renderFrameTimed(cam, outHeight, outWidth, frame::colormap::RGB, invFrameRate, reflectCount, globalIllumination, useLod);
         } else {
-            currentPreviewFrame = sim.grid.renderFrameTimed(cam, outHeight, outWidth, frame::colormap::RGB);
+            currentPreviewFrame = sim.grid.fastRenderFrame(cam, outHeight, outWidth, frame::colormap::RGB);
         }
         
         if (textu == 0) {

util/grid/grid3eigen.hpp

@@ -21,6 +21,7 @@
 #include <unordered_set>
 #include <random>
 #include <chrono>
+#include <cstdint>
 
 #ifdef SSE
 #include <immintrin.h>
@@ -28,7 +29,7 @@
 
 constexpr int Dim = 3;
 
-template<typename T>
+template<typename T, typename IndexType = uint16_t>
 class Octree {
 public:
     using PointType = Eigen::Matrix<float, Dim, 1>;
@@ -40,37 +41,48 @@ public:
         }
     };
 
-    struct NodeData {
+    struct Material {
-        T data;
-        PointType position;
-        int objectId;
-        int subId;
-        bool active;
-        bool visible;
-        float size;
-        Eigen::Vector3f color;
         float emittance;
         float roughness;
         float metallic;
         float transmission;
         float ior;
 
-        NodeData(const T& data, const PointType& pos, bool visible, Eigen::Vector3f color, float size = 0.01f,
+        Material(float e = 0.0f, float r = 1.0f, float m = 0.0f, float t = 0.0f, float i = 1.45f)
-                 bool active = true, int objectId = -1, int subId = 0, float emittance = 0.0f, 
+            : emittance(e), roughness(r), metallic(m), transmission(t), ior(i) {}
-                 float roughness = 1.0f, float metallic = 0.0f, float transmission = 0.0f, float ior = 1.45f) 
-                : data(data), position(pos), objectId(objectId), subId(subId), active(active), visible(visible), 
-                  color(color), size(size), emittance(emittance), roughness(roughness), metallic(metallic), 
-                  transmission(transmission), ior(ior) {}
 
-        NodeData() : objectId(-1), subId(0), active(false), visible(false), size(0.0f), emittance(0.0f), 
+        bool operator<(const Material& o) const {
-                     roughness(1.0f), metallic(0.0f), transmission(0.0f), ior(1.45f) {}
+            if (emittance != o.emittance) return emittance < o.emittance;
+            if (roughness != o.roughness) return roughness < o.roughness;
+            if (metallic != o.metallic) return metallic < o.metallic;
+            if (transmission != o.transmission) return transmission < o.transmission;
+            return ior < o.ior;
+        }
+    };
+    
+    struct NodeData {
+        T data;
+        PointType position;
+        int objectId;
+        int subId;
+        float size;
+        IndexType colorIdx;
+        IndexType materialIdx;
+        bool active;
+        bool visible;
+
+        NodeData(const T& data, const PointType& pos, bool visible, IndexType colorIdx, float size = 0.01f,
+                 bool active = true, int objectId = -1, int subId = 0, IndexType materialIdx = 0) 
+                : data(data), position(pos), objectId(objectId), subId(subId), size(size), 
+                  colorIdx(colorIdx), materialIdx(materialIdx), active(active), visible(visible) {}
+        
+        NodeData() : objectId(-1), subId(0), size(0.0f), colorIdx(0), materialIdx(0), 
+                     active(false), visible(false) {}
         
-        // Helper method to get half-size for cube
         PointType getHalfSize() const {
             return PointType(size * 0.5f, size * 0.5f, size * 0.5f);
         }
         
-        // Helper method to get bounding box for cube
         BoundingBox getCubeBounds() const {
             PointType halfSize = getHalfSize();
             return {position - halfSize, position + halfSize};
@@ -116,10 +128,85 @@ private:
     Eigen::Vector3f skylight_ = {0.1f, 0.1f, 0.1f};
     Eigen::Vector3f backgroundColor_ = {0.53f, 0.81f, 0.92f};
     
+    // Addressable Maps
+    std::unique_ptr<std::mutex> mapMutex_;
+    std::vector<Eigen::Vector3f> colorMap_;
+    std::map<Eigen::Vector3f, IndexType, Vector3fCompare> colorToIndex_;
+    
+    std::vector<Material> materialMap_;
+    std::map<Material, IndexType> materialToIndex_;
+
     std::map<std::pair<int, int>, std::shared_ptr<Mesh>> meshCache_;
     std::set<std::pair<int, int>> dirtyMeshes_;
     int nextSubIdGenerator_ = 1;
     
+public:
+    inline IndexType getColorIndex(const Eigen::Vector3f& color) {
+        std::lock_guard<std::mutex> lock(*mapMutex_);
+        auto it = colorToIndex_.find(color);
+        if (it != colorToIndex_.end()) return it->second;
+        
+        if (colorMap_.size() >= std::numeric_limits<IndexType>::max()) {
+            IndexType bestIdx = 0;
+            float bestDist = std::numeric_limits<float>::max();
+            for (size_t i = 0; i < colorMap_.size(); ++i) {
+                float dist = (colorMap_[i] - color).squaredNorm();
+                if (dist < bestDist) {
+                    bestDist = dist;
+                    bestIdx = static_cast<IndexType>(i);
+                }
+            }
+            return bestIdx;
+        }
+
+        IndexType idx = static_cast<IndexType>(colorMap_.size());
+        colorMap_.push_back(color);
+        colorToIndex_[color] = idx;
+        return idx;
+    }
+
+    inline const Eigen::Vector3f& getColor(IndexType idx) const {
+        if (idx < colorMap_.size()) return colorMap_[idx];
+        static const Eigen::Vector3f fallback = Eigen::Vector3f::Zero();
+        return fallback;
+    }
+
+    inline IndexType getMaterialIndex(const Material& mat) {
+        std::lock_guard<std::mutex> lock(*mapMutex_);
+        auto it = materialToIndex_.find(mat);
+        if (it != materialToIndex_.end()) return it->second;
+
+        if (materialMap_.size() >= std::numeric_limits<IndexType>::max()) {
+            IndexType bestIdx = 0;
+            float bestDist = std::numeric_limits<float>::max();
+            for (size_t i = 0; i < materialMap_.size(); ++i) {
+                float d_e = materialMap_[i].emittance - mat.emittance;
+                float d_r = materialMap_[i].roughness - mat.roughness;
+                float d_m = materialMap_[i].metallic - mat.metallic;
+                float d_t = materialMap_[i].transmission - mat.transmission;
+                float d_i = materialMap_[i].ior - mat.ior;
+                float dist = d_e*d_e + d_r*d_r + d_m*d_m + d_t*d_t + d_i*d_i;
+                if (dist < bestDist) {
+                    bestDist = dist;
+                    bestIdx = static_cast<IndexType>(i);
+                }
+            }
+            return bestIdx;
+        }
+
+        IndexType idx = static_cast<IndexType>(materialMap_.size());
+        materialMap_.push_back(mat);
+        materialToIndex_[mat] = idx;
+        return idx;
+    }
+
+    inline const Material& getMaterial(IndexType idx) const {
+        if (idx < materialMap_.size()) return materialMap_[idx];
+        static const Material fallback;
+        return fallback;
+    }
+
+private:
     void invalidateMesh(int objectId, int subId) {
         if (objectId < 0) return;
         dirtyMeshes_.insert({objectId, subId});
@@ -272,7 +359,7 @@ private:
         }
     }
     
-    void ensureLOD(const OctreeNode* node) const {
+    void ensureLOD(OctreeNode* node) {
         std::lock_guard<std::mutex> lock(node->lodMutex);
         if (node->lodData != nullptr) return;
 
@@ -294,12 +381,13 @@ private:
 
         auto accumulate = [&](const std::shared_ptr<NodeData>& item) {
             if (!item || !item->active || !item->visible) return;
-            avgColor += item->color;
+            avgColor += getColor(item->colorIdx);
-            avgEmittance += item->emittance;
+            Material mat = getMaterial(item->materialIdx);
-            avgRoughness += item->roughness;
+            avgEmittance += mat.emittance;
-            avgMetallic += item->metallic;
+            avgRoughness += mat.roughness;
-            avgTransmission += item->transmission;
+            avgMetallic += mat.metallic;
-            avgIor += item->ior;
+            avgTransmission += mat.transmission;
+            avgIor += mat.ior;
             count++;
         };
 
@@ -318,18 +406,16 @@ private:
             float invCount = 1.0f / count;
             
             auto lod = std::make_shared<NodeData>();
-            
             lod->position = node->center;
-            lod->color = avgColor * invCount;
             
             PointType nodeDims = node->bounds.second - node->bounds.first;
             lod->size = nodeDims.maxCoeff();
 
-            lod->emittance = avgEmittance * invCount;
+            lod->colorIdx = getColorIndex(avgColor * invCount);
-            lod->roughness = avgRoughness * invCount;
+            Material avgMat(avgEmittance * invCount, avgRoughness * invCount, 
-            lod->metallic = avgMetallic * invCount;
+                            avgMetallic * invCount, avgTransmission * invCount, avgIor * invCount);
-            lod->transmission = avgTransmission * invCount;
+            lod->materialIdx = getMaterialIndex(avgMat);
-            lod->ior = avgIor * invCount;
+            
             lod->active = true;
             lod->visible = true;
             lod->objectId = -1; 
@@ -546,14 +632,17 @@ private:
         Ray ray(rayOrig, rayDir);
         rayCubeIntersect(ray, obj.get(), t, normal, hitPoint);
         
+        Eigen::Vector3f objColor = getColor(obj->colorIdx);
+        Material objMat = getMaterial(obj->materialIdx);
+
         Eigen::Vector3f finalColor = globalIllumination ? skylight_ : Eigen::Vector3f::Zero();
-        if (obj->emittance > 0.0f) {
+        if (objMat.emittance > 0.0f) {
-            finalColor += obj->color * obj->emittance;
+            finalColor += objColor * objMat.emittance;
         }
 
-        float roughness = std::clamp(obj->roughness, 0.01f, 1.0f);
+        float roughness = std::clamp(objMat.roughness, 0.01f, 1.0f);
-        float metallic = std::clamp(obj->metallic, 0.0f, 1.0f);
+        float metallic = std::clamp(objMat.metallic, 0.0f, 1.0f);
-        float transmission = std::clamp(obj->transmission, 0.0f, 1.0f);
+        float transmission = std::clamp(objMat.transmission, 0.0f, 1.0f);
         
         PointType V = -rayDir;
         float cosThetaI = normal.dot(V);
@@ -565,7 +654,7 @@ private:
         float rayOffset = std::max(1e-4f, 1e-5f * coordMax);
 
         Eigen::Vector3f F0 = Eigen::Vector3f::Constant(0.04f);
-        F0 = F0 * (1.0f - metallic) + obj->color * metallic;
+        F0 = F0 * (1.0f - metallic) + objColor * metallic;
 
         PointType H = sampleGGX(n_eff, roughness, rngState);
         float VdotH = std::max(0.001f, V.dot(H));
@@ -594,23 +683,23 @@ private:
         float diffuseWeight = (1.0f - transmission) * (1.0f - metallic);
 
         if (transmissionWeight > 0.0f) {
-            float eta = isInside ? obj->ior : (1.0f / obj->ior);
+            float eta = isInside ? objMat.ior : (1.0f / objMat.ior);
             float k = 1.0f - eta * eta * (1.0f - VdotH * VdotH);
             
             if (k >= 0.0f) {
                 secondDir = ((eta * VdotH - std::sqrt(k)) * H - eta * V).normalized();
                 secondOrigin = hitPoint - n_eff * rayOffset;
                 W_second = (Eigen::Vector3f::Constant(1.0f) - F_spec) * transmissionWeight;
-                W_second = W_second.cwiseProduct(obj->color);
+                W_second = W_second.cwiseProduct(objColor);
             } else {
                 Eigen::Vector3f tirWeight = (Eigen::Vector3f::Constant(1.0f) - F_spec) * transmissionWeight;
-                W_spec += tirWeight.cwiseProduct(obj->color);
+                W_spec += tirWeight.cwiseProduct(objColor);
             }
         } else if (diffuseWeight > 0.0f) {
             secondDir = sampleCosineHemisphere(n_eff, rngState);
             secondOrigin = hitPoint + n_eff * rayOffset;
             W_second = (Eigen::Vector3f::Constant(1.0f) - F_spec) * diffuseWeight;
-            W_second = W_second.cwiseProduct(obj->color);
+            W_second = W_second.cwiseProduct(objColor);
         }
 
         W_spec = W_spec.cwiseMin(Eigen::Vector3f::Constant(4.0f));
@@ -738,22 +827,22 @@ private:
     }
 
     template<typename V>
-    void writeVal(std::ofstream& out, const V& val) const {
+    inline void writeVal(std::ofstream& out, const V& val) const {
         out.write(reinterpret_cast<const char*>(&val), sizeof(V));
     }
 
     template<typename V>
-    void readVal(std::ifstream& in, V& val) {
+    inline void readVal(std::ifstream& in, V& val) {
         in.read(reinterpret_cast<char*>(&val), sizeof(V));
     }
 
-    void writeVec3(std::ofstream& out, const Eigen::Vector3f& vec) const {
+    inline void writeVec3(std::ofstream& out, const Eigen::Vector3f& vec) const {
         writeVal(out, vec.x());
         writeVal(out, vec.y());
         writeVal(out, vec.z());
     }
 
-    void readVec3(std::ifstream& in, Eigen::Vector3f& vec) {
+    inline void readVec3(std::ifstream& in, Eigen::Vector3f& vec) {
         float x, y, z;
         readVal(in, x); readVal(in, y); readVal(in, z);
         vec = Eigen::Vector3f(x, y, z);
@@ -766,20 +855,14 @@ private:
             size_t pointCount = node->points.size();
             writeVal(out, pointCount);
             for (const auto& pt : node->points) {
-                // Write raw data T (Must be POD)
                 writeVal(out, pt->data);
-                // Write properties
                 writeVec3(out, pt->position);
                 writeVal(out, pt->objectId);
                 writeVal(out, pt->active);
                 writeVal(out, pt->visible);
                 writeVal(out, pt->size);
-                writeVec3(out, pt->color);
+                writeVal(out, pt->colorIdx);
-                writeVal(out, pt->emittance);
+                writeVal(out, pt->materialIdx);
-                writeVal(out, pt->roughness);
-                writeVal(out, pt->metallic);
-                writeVal(out, pt->transmission);
-                writeVal(out, pt->ior);
             }
         } else {
             // Write bitmask of active children
@@ -818,12 +901,8 @@ private:
             readVal(in, pt->active);
             readVal(in, pt->visible);
             readVal(in, pt->size);
-            readVec3(in, pt->color);
+            readVal(in, pt->colorIdx);
-            readVal(in, pt->emittance);
+            readVal(in, pt->materialIdx);
-            readVal(in, pt->roughness);
-            readVal(in, pt->metallic);
-            readVal(in, pt->transmission);
-            readVal(in, pt->ior);
             node->points.push_back(pt);
         }
 
@@ -1093,7 +1172,7 @@ private:
                 float x = distSq / rSq;
                 float w = (1.0f - x) * (1.0f - x);
                 density += w;
-                accumulatedColor += neighbor->color * w;
+                accumulatedColor += getColor(neighbor->colorIdx) * w;
                 totalWeight += w;
             }
         }
@@ -1145,9 +1224,9 @@ private:
 public:
     Octree(const PointType& minBound, const PointType& maxBound, size_t maxPointsPerNode=8, size_t maxDepth = 16) :
     root_(std::make_unique<OctreeNode>(minBound, maxBound)), maxPointsPerNode(maxPointsPerNode),
-    maxDepth(maxDepth), size(0) {}
+    maxDepth(maxDepth), size(0), mapMutex_(std::make_unique<std::mutex>()) {}
 
-    Octree() : root_(nullptr), maxPointsPerNode(8), maxDepth(16), size(0) {}
+    Octree() : root_(nullptr), maxPointsPerNode(8), maxDepth(16), size(0), mapMutex_(std::make_unique<std::mutex>()) {}
 
     void setSkylight(const Eigen::Vector3f& skylight) { 
         skylight_ = skylight; 
@@ -1177,8 +1256,13 @@ public:
     bool set(const T& data, const PointType& pos, bool visible, Eigen::Vector3f color, float size = 0.01f, bool active = true,
              int objectId = -1, int subId = 0, float emittance = 0.0f, float roughness = 1.0f, 
              float metallic = 0.0f, float transmission = 0.0f, float ior = 1.45f) {
-        auto pointData = std::make_shared<NodeData>(data, pos, visible, color, size, active, objectId, subId,
+        
-                                                    emittance, roughness, metallic, transmission, ior);
+        IndexType cIdx = getColorIndex(color);
+        Material mat(emittance, roughness, metallic, transmission, ior);
+        IndexType mIdx = getMaterialIndex(mat);
+        
+        auto pointData = std::make_shared<NodeData>(data, pos, visible, cIdx, size, active, objectId, subId, mIdx);
+        
         if (insertRecursive(root_.get(), pointData, 0)) {
             return true;
         }
@@ -1196,10 +1280,28 @@ public:
         writeVal(out, maxPointsPerNode);
         writeVal(out, size);
         
-        // Save global settings
         writeVec3(out, skylight_);
         writeVec3(out, backgroundColor_);
 
+        {
+            std::lock_guard<std::mutex> lock(*mapMutex_);
+            size_t cMapSize = colorMap_.size();
+            writeVal(out, cMapSize);
+            for (const auto& c : colorMap_) {
+                writeVec3(out, c);
+            }
+
+            size_t mMapSize = materialMap_.size();
+            writeVal(out, mMapSize);
+            for (const auto& m : materialMap_) {
+                writeVal(out, m.emittance);
+                writeVal(out, m.roughness);
+                writeVal(out, m.metallic);
+                writeVal(out, m.transmission);
+                writeVal(out, m.ior);
+            }
+        }
+        
         writeVec3(out, root_->bounds.first);
         writeVec3(out, root_->bounds.second);
 
@@ -1228,13 +1330,40 @@ public:
         readVec3(in, skylight_);
         readVec3(in, backgroundColor_);
 
+        {
+            std::lock_guard<std::mutex> lock(*mapMutex_);
+            colorMap_.clear();
+            colorToIndex_.clear();
+            materialMap_.clear();
+            materialToIndex_.clear();
+
+            size_t cMapSize;
+            readVal(in, cMapSize);
+            colorMap_.resize(cMapSize);
+            for (size_t i = 0; i < cMapSize; ++i) {
+                readVec3(in, colorMap_[i]);
+                colorToIndex_[colorMap_[i]] = static_cast<IndexType>(i);
+            }
+
+            size_t mMapSize;
+            readVal(in, mMapSize);
+            materialMap_.resize(mMapSize);
+            for (size_t i = 0; i < mMapSize; ++i) {
+                readVal(in, materialMap_[i].emittance);
+                readVal(in, materialMap_[i].roughness);
+                readVal(in, materialMap_[i].metallic);
+                readVal(in, materialMap_[i].transmission);
+                readVal(in, materialMap_[i].ior);
+                materialToIndex_[materialMap_[i]] = static_cast<IndexType>(i);
+            }
+        }
+
         PointType minBound, maxBound;
         readVec3(in, minBound);
         readVec3(in, maxBound);
 
         root_ = std::make_unique<OctreeNode>(minBound, maxBound);
-        std::multimap<Eigen::Vector3f, std::shared_ptr<NodeData>, Vector3fCompare> pointMap;
+        deserializeNode(in, root_.get());
-        deserializeNode(in, root_.get(), pointMap);
 
         in.close();
         std::cout << "successfully loaded grid from " << filename << std::endl;
@@ -1300,17 +1429,36 @@ public:
         pointData->data = newData;
         pointData->position = newPos;
         pointData->visible = newVisible;
-        if (newColor != Eigen::Vector3f(1.0f, 1.0f, 1.0f)) pointData->color = newColor;
+        
+        if (newColor != Eigen::Vector3f(1.0f, 1.0f, 1.0f)) pointData->colorIdx = getColorIndex(newColor);
         if (newSize > 0) pointData->size = newSize;
         pointData->active = newActive;
         pointData->objectId = targetObjId;
         pointData->subId = finalSubId;
         
-        if (newEmittance >= 0) pointData->emittance = newEmittance;
+        Material mat = getMaterial(pointData->materialIdx);
-        if (newRoughness >= 0) pointData->roughness = newRoughness;
+        bool matChanged = false;
-        if (newMetallic >= 0) pointData->metallic = newMetallic;
+        if (newEmittance >= 0) {
-        if (newTransmission >= 0) pointData->transmission = newTransmission;
+            mat.emittance = newEmittance;
-        if (newIor >= 0) pointData->ior = newIor;
+            matChanged = true;
+        }
+        if (newRoughness >= 0) { 
+            mat.roughness = newRoughness;
+            matChanged = true;
+        }
+        if (newMetallic >= 0) { 
+            mat.metallic = newMetallic;
+            matChanged = true;
+        }
+        if (newTransmission >= 0) { 
+            mat.transmission = newTransmission;
+            matChanged = true;
+        }
+        if (newIor >= 0) { 
+            mat.ior = newIor;
+            matChanged = true;
+        }
+        if (matChanged) pointData->materialIdx = getMaterialIndex(mat);
         
         bool res = insertRecursive(root_.get(), pointData, 0);
         
@@ -1326,9 +1474,7 @@ public:
         auto pointData = find(pos);
         if (!pointData) return false;
 
-        bool sizeDecremented = false;
         removeRecursive(root_.get(), pointData->getCubeBounds(), pos, EPSILON);
-        
         pointData->position = newPos;
 
         if (insertRecursive(root_.get(), pointData, 0)) {
@@ -1372,7 +1518,7 @@ public:
     bool setColor(const PointType& pos, Eigen::Vector3f color, float tolerance = EPSILON) {
         auto pointData = find(pos, tolerance);
         if (!pointData) return false;
-        pointData->color = color;
+        pointData->colorIdx = getColorIndex(color);
         invalidateLODForPoint(pointData);
         return true;
     }
@@ -1380,7 +1526,9 @@ public:
     bool setEmittance(const PointType& pos, float emittance, float tolerance = EPSILON) {
         auto pointData = find(pos, tolerance);
         if (!pointData) return false;
-        pointData->emittance = emittance;
+        Material mat = getMaterial(pointData->materialIdx);
+        mat.emittance = emittance;
+        pointData->materialIdx = getMaterialIndex(mat);
         invalidateLODForPoint(pointData);
         return true;
     }
@@ -1388,7 +1536,9 @@ public:
     bool setRoughness(const PointType& pos, float roughness, float tolerance = EPSILON) {
         auto pointData = find(pos, tolerance);
         if (!pointData) return false;
-        pointData->roughness = roughness;
+        Material mat = getMaterial(pointData->materialIdx);
+        mat.roughness = roughness;
+        pointData->materialIdx = getMaterialIndex(mat);
         invalidateLODForPoint(pointData);
         return true;
     }
@@ -1396,7 +1546,9 @@ public:
     bool setMetallic(const PointType& pos, float metallic, float tolerance = EPSILON) {
         auto pointData = find(pos, tolerance);
         if (!pointData) return false;
-        pointData->metallic = metallic;
+        Material mat = getMaterial(pointData->materialIdx);
+        mat.metallic = metallic;
+        pointData->materialIdx = getMaterialIndex(mat);
         invalidateLODForPoint(pointData);
         return true;
     }
@@ -1404,7 +1556,9 @@ public:
     bool setTransmission(const PointType& pos, float transmission, float tolerance = EPSILON) {
         auto pointData = find(pos, tolerance);
         if (!pointData) return false;
-        pointData->transmission = transmission;
+        Material mat = getMaterial(pointData->materialIdx);
+        mat.transmission = transmission;
+        pointData->materialIdx = getMaterialIndex(mat);
         invalidateLODForPoint(pointData);
         return true;
     }
@@ -1533,10 +1687,11 @@ public:
                     PointType normal, hitPoint;
 
                     rayCubeIntersect(ray, obj.get(), t, normal, hitPoint);
-                    color = obj->color;
+                    color = getColor(obj->colorIdx);
+                    Material objMat = getMaterial(obj->materialIdx);
                     
-                    if (obj->emittance > 0.0f) {
+                    if (objMat.emittance > 0.0f) {
-                        color = color * obj->emittance;
+                        color = color * objMat.emittance;
                     } else {
                         float diffuse = std::max(0.0f, normal.dot(globalLightDir));
                         float ambient = 0.35f;
@@ -1607,10 +1762,11 @@ public:
                     PointType normal, hitPoint;
 
                     rayCubeIntersect(ray, obj.get(), t, normal, hitPoint);
-                    color = obj->color;
+                    color = getColor(obj->colorIdx);
+                    Material objMat = getMaterial(obj->materialIdx);
                     
-                    if (obj->emittance > 0.0f) {
+                    if (objMat.emittance > 0.0f) {
-                        color = color * obj->emittance;
+                        color = color * objMat.emittance;
                     } else {
                         float diffuse = std::max(0.0f, normal.dot(globalLightDir));
                         float ambient = 0.35f;
@@ -1847,7 +2003,6 @@ public:
         std::vector<Eigen::Vector3f> vertexColors;
         std::vector<Eigen::Vector3i> triangles;
 
-        // Marching Cubes Loop
         for(int z = 0; z < resolution; ++z) {
             for(int y = 0; y < resolution; ++y) {
                 for(int x = 0; x < resolution; ++x) {
@@ -1918,10 +2073,13 @@ public:
             PointType dirs[6] = {{1,0,0}, {-1,0,0}, {0,1,0}, {0,-1,0}, {0,0,1}, {0,0,-1}};
             for(int i=0; i<6; ++i) {
                 auto neighbor = find(node->position + dirs[i] * node->size, checkRad);
-                if(neighbor && neighbor->objectId == objectId && neighbor->active && neighbor->transmission < 0.01f) {
+                if(neighbor && neighbor->objectId == objectId && neighbor->active) {
+                    Material nMat = getMaterial(neighbor->materialIdx);
+                    if (nMat.transmission < 0.01f) {
                         hiddenSides++;
                     }
                 }
+            }
 
             if(hiddenSides == 6) {
                 if(node->subId != -999) {
@@ -1981,6 +2139,8 @@ public:
         
         size_t nodeMem = totalNodes * sizeof(OctreeNode);
         size_t dataMem = actualPoints * (sizeof(NodeData) + 16); 
+        size_t mapMem = colorMap_.size() * sizeof(Eigen::Vector3f) + materialMap_.size() * sizeof(Material);
+        size_t maxSize = ((1 << (sizeof(IndexType)*8 - 2) - 1) * 2) + 1;
 
         os << "========================================\n";
         os << "             OCTREE STATS               \n";
@@ -2001,12 +2161,16 @@ public:
         os << "  Points/Leaf (Avg) : " << std::fixed << std::setprecision(2) << avgPointsPerLeaf << "\n";
         os << "  Points/Leaf (Min) : " << minPointsInLeaf << "\n";
         os << "  Points/Leaf (Max) : " << maxPointsInLeaf << "\n";
+        os << "Maps:\n";
+        os << "  Unique Colors     : " << colorMap_.size() << "/" << maxSize << "\n";
+        os << "  Unique Materials  : " << materialMap_.size() << "/" << maxSize << "\n";
         os << "Bounds:\n";
         os << "  Min               : [" << root_->bounds.first.transpose() << "]\n";
         os << "  Max               : [" << root_->bounds.second.transpose() << "]\n";
         os << "Memory (Approx):\n";
         os << "  Node Structure    : " << (nodeMem / 1024.0) << " KB\n";
         os << "  Point Data        : " << (dataMem / 1024.0) << " KB\n";
+        os << "  Dictionary Maps   : " << (mapMem / 1024.0) << " KB\n";
         os << "========================================\n" << std::defaultfloat;
     }
 
@@ -2020,6 +2184,14 @@ public:
         PointType maxBound = root_->bounds.second;
         root_ = std::make_unique<OctreeNode>(minBound, maxBound);
         
+        {
+            std::lock_guard<std::mutex> lock(*mapMutex_);
+            colorMap_.clear();
+            colorToIndex_.clear();
+            materialMap_.clear();
+            materialToIndex_.clear();
+        }
+        
         size = 0;
     }
 };

util/sim/planet.hpp

@@ -13,6 +13,7 @@
 #include <algorithm>
 #include <queue>
 #include <unordered_map>
+#include <set>
 
 #include "../grid/grid3eigen.hpp"
 #include "../timing_decorator.cpp"
@@ -522,12 +523,10 @@ public:
     void extraplateste() {
         TIME_FUNCTION;
         std::uniform_real_distribution<float> distFloat(0.0f, 1.0f);
-        
+        std::vector<std::pair<int, float>> plateStats;
-        struct PlateStats {
+        for (int i = 0; i < config.numPlates; i++) {
-            int id;
+            plateStats.emplace_back(std::make_pair(i, 0.0f));
-            float avgElevation;
+        }
-        };
-        std::vector<PlateStats> stats(config.numPlates);
 
         for (int i = 0; i < config.numPlates; i++) {
             float sumElevation = 0.0f;
@@ -539,7 +538,7 @@ public:
             }
             
             if (!plates[i].assignedNodes.empty()) {
-                stats[i].avgElevation = sumElevation / plates[i].assignedNodes.size();
+                plateStats[i].first = sumElevation / plates[i].assignedNodes.size();
                 centroid /= plates[i].assignedNodes.size();
 
                 float maxSpread = 0.0f;
@@ -564,9 +563,8 @@ public:
                     plates[i].plateEulerPole = config.surfaceNodes[bestNodeIdx]; 
                 }
             } else {
-                stats[i].avgElevation = config.radius;
+                plateStats[i].first = config.radius;
             }
-            stats[i].id = i;
             
             Eigen::Vector3f randomDir(distFloat(rng) - 0.5f, distFloat(rng) - 0.5f, distFloat(rng) - 0.5f);
             randomDir.normalize();
@@ -579,15 +577,15 @@ public:
             plates[i].temperature = distFloat(rng) * 1000.0f;
         }
         
-        std::sort(stats.begin(), stats.end(), [](const PlateStats& a, const PlateStats& b) {
+        std::sort(plateStats.begin(), plateStats.end(), [](const std::pair<int, float>& a, const std::pair<int, float>& b) {
-            return a.avgElevation < b.avgElevation;
+            return a.second < b.second;
         });
         
         int oneThird = config.numPlates / 3;
         int twoThirds = (2 * config.numPlates) / 3;
         
         for (int i = 0; i < config.numPlates; i++) {
-            int pID = stats[i].id;
+            int pID = plateStats[i].first;
             if (i < oneThird) {
                 plates[pID].ptype = PlateType::OCEANIC;
                 plates[pID].thickness = distFloat(rng) * 10.0f + 5.0f;
@@ -773,7 +771,6 @@ public:
                     
                     if (w1 > 0.99f || w2 > 0.99f || w3 > 0.99f) continue;
 
-                    // Calculate position
                     v3 interpNormal = (p1.originalPos.cast<float>().normalized() * w1 + p2.originalPos.cast<float>().normalized() * w2 + p3.originalPos.cast<float>().normalized() * w3);
                     interpNormal.normalize(); 
                     
@@ -792,7 +789,6 @@ public:
                     newPt.surface = true;
                     newPt.currentPos = smoothPos;
                     
-                    // Assign properties based on dominant weight
                     if (w1 > w2 && w1 > w3) {
                         newPt.plateID = p1.plateID;
                         newPt.originColor = p1.originColor;
@@ -808,10 +804,8 @@ public:
                     newPt.noisePos = (p1.noisePos.cast<float>() * w1 + p2.noisePos.cast<float>() * w2 + p3.noisePos.cast<float>() * w3).cast<Eigen::half>();
                     newPt.tectonicPos = (p1.tectonicPos.cast<float>() * w1 + p2.tectonicPos.cast<float>() * w2 + p3.tectonicPos.cast<float>() * w3).cast<Eigen::half>();
 
-                    // Insert into Grid
                     grid.set(newPt, newPt.currentPos, true, newPt.originColor.cast<float>(), config.voxelSize, true, 1, 2, false, 0.0f, 0.0f, 0.0f);
                     
-                    // FIX: Save to interpolatedNodes so we don't lose the data reference
                     config.interpolatedNodes.push_back(newPt);
 
                     counter++;
@@ -820,233 +814,60 @@ public:
         }
         grid.optimize();
         std::cout << "Interpolated " << counter << " surface gaps." << std::endl;
-        
-        sealCracks();
-    }
-
-    void sealCracks() {
-        TIME_FUNCTION;
-        std::vector<Particle> patchNodes;
-        float vsize = config.voxelSize;
-
-        std::vector<Particle*> candidates;
-        candidates.reserve(config.interpolatedNodes.size() + config.surfaceNodes.size());
-        for(auto& p : config.surfaceNodes) candidates.push_back(&p);
-        for(auto& p : config.interpolatedNodes) candidates.push_back(&p);
-
-        int patchedCount = 0;
-        
-        std::vector<v3> directions = {
-            v3(vsize,0,0), v3(-vsize,0,0),
-            v3(0,vsize,0), v3(0,-vsize,0),
-            v3(0,0,vsize), v3(0,0,-vsize)
-        };
-
-        for (Particle* p : candidates) {
-            for (const auto& dir : directions) {
-                v3 checkPos = p->currentPos + dir;
-                
-                if (grid.find(checkPos, vsize * 0.1f) == nullptr) {
-                    
-                    int neighborsFound = 0;
-                    for (const auto& nDir : directions) {
-                        if (grid.find(checkPos + nDir, vsize * 0.1f) != nullptr) {
-                            neighborsFound++;
-                        }
-                    }
-
-                    if (neighborsFound >= 4) {
-                        Particle patch = *p;
-                        patch.currentPos = checkPos;
-                        patch.tectonicPos = checkPos.cast<Eigen::half>();
-                        
-                        bool alreadyPatched = false;
-                        for(const auto& pp : patchNodes) {
-                            if((pp.currentPos - checkPos).norm() < 1.0f) { alreadyPatched=true; break; }
-                        }
-                        
-                        if (!alreadyPatched) {
-                            patchNodes.push_back(patch);
-                            grid.set(patch, checkPos, true, patch.originColor.cast<float>(), vsize, true, 1, 2, false, 0.0f, 0.0f, 0.0f);
-                            patchedCount++;
-                        }
-                    }
-                }
-            }
-        }
-
-        for(const auto& p : patchNodes) {
-            config.interpolatedNodes.push_back(p);
-        }
-        
-        std::cout << "Sealed " << patchedCount << " surface cracks." << std::endl;
-        grid.optimize();
     }
 
     void fillPlanet() {
         TIME_FUNCTION;
-        std::cout << "Starting Volume Fill (Naive)..." << std::endl;
+        if (config.interpolatedNodes.empty()) {
-
+            std::cout << "Please run interpolate surface first." << std::endl;
-        std::vector<Particle*> hullPtrs;
+            return;
-        hullPtrs.reserve(config.surfaceNodes.size() + config.interpolatedNodes.size());
-        
-        float maxDistSq = 0.0f;
-        for (size_t i = 0; i < config.surfaceNodes.size(); i++) {
-            hullPtrs.push_back(&config.surfaceNodes[i]);
-            float d2 = config.surfaceNodes[i].currentPos.squaredNorm();
-            if (d2 > maxDistSq) maxDistSq = d2;
         }
-        for (size_t i = 0; i < config.interpolatedNodes.size(); i++) {
+        std::cout << "Starting Volume Fill..." << std::endl;
-            hullPtrs.push_back(&config.interpolatedNodes[i]);
+
-            float d2 = config.interpolatedNodes[i].currentPos.squaredNorm();
+        float safeRadius = config.radius - std::abs(config.valleyDepth) - (config.noiseStrength * 2.0f) - config.voxelSize;
-            if (d2 > maxDistSq) maxDistSq = d2;
+        if (safeRadius <= 0) safeRadius = config.radius * 0.5f;
+
+        int maxSteps = std::ceil(safeRadius / config.voxelSize);
+        size_t fillCount = 0;
+
+        for (int x = -maxSteps; x <= maxSteps; ++x) {
+            for (int y = -maxSteps; y <= maxSteps; ++y) {
+                for (int z = -maxSteps; z <= maxSteps; ++z) {
+                    v3 pos = config.center + v3(x, y, z) * config.voxelSize;
+                    float dist = (pos - config.center).norm();
+
+                    if (dist <= safeRadius) {
+                        if (grid.find(pos, config.voxelSize * 0.5f) == nullptr) {
+                            Particle ip;
+                            ip.surface = false;
+                            ip.plateID = -1;
+                            ip.currentPos = pos;
+                            ip.originalPos = pos.cast<Eigen::half>();
+                            ip.noisePos = pos.cast<Eigen::half>();
+                            ip.tectonicPos = pos.cast<Eigen::half>();
+
+                            float depthRatio = dist / safeRadius;
+                            Eigen::Vector3f coreColor(1.0f, 0.9f, 0.4f);
+                            Eigen::Vector3f mantleColor(0.8f, 0.15f, 0.0f);
+                            Eigen::Vector3f finalColor = mantleColor;
+
+                            if (depthRatio < 0.5f) {
+                                float blend = depthRatio * 2.0f;
+                                finalColor = coreColor * (1.0f - blend) + mantleColor * blend;
                             }
 
-        float maxRadius = std::sqrt(maxDistSq);
+                            ip.originColor = finalColor.cast<Eigen::half>();
-        float step = config.voxelSize * 0.5f; 
 
-        float cellScale = 1.0f / (config.voxelSize * 2.0f);
+                            grid.set(ip, pos, true, finalColor, config.voxelSize, true, 1, 3, false, 0.0f, 0.0f, 0.0f);
-        int tableSize = hullPtrs.size() * 2 + 1;
+                            fillCount++;
-        
-        std::vector<int> head(tableSize, -1);
-        std::vector<int> next(hullPtrs.size(), -1);
-
-        for (int i = 0; i < hullPtrs.size(); i++) {
-            v3 p = hullPtrs[i]->currentPos;
-            int64_t x = static_cast<int64_t>(std::floor(p.x() * cellScale));
-            int64_t y = static_cast<int64_t>(std::floor(p.y() * cellScale));
-            int64_t z = static_cast<int64_t>(std::floor(p.z() * cellScale));
-            
-            uint64_t h = (x * 73856093) ^ (y * 19349663) ^ (z * 83492791);
-            int bucket = h % tableSize;
-            
-            next[i] = head[bucket];
-            head[bucket] = i;
-        }
-
-        std::cout << "Spatial Map Built. Scanning volume..." << std::endl;
-
-        struct ThreadData {
-            std::vector<Particle> gridCandidates;
-            std::vector<Particle> surfaceCandidates;
-        };
-        
-        int gridLimit = static_cast<int>(std::ceil(maxRadius / step));
-
-        #pragma omp parallel
-        {
-            ThreadData tData;
-            
-            #pragma omp for collapse(3) schedule(dynamic, 8)
-            for (int x = -gridLimit; x <= gridLimit; x++) {
-                for (int y = -gridLimit; y <= gridLimit; y++) {
-                    for (int z = -gridLimit; z <= gridLimit; z++) {
-                        
-                        v3 pos(x * step, y * step, z * step);
-                        float dCentSq = pos.squaredNorm();
-
-                        if (dCentSq > maxDistSq) continue;
-
-                        if (grid.find(pos, step * 0.1f) != nullptr) continue;
-
-                        Particle* nearest[3] = {nullptr, nullptr, nullptr};
-                        float dists[3] = {1e20f, 1e20f, 1e20f};
-                        int foundCount = 0;
-
-                        int64_t bx = static_cast<int64_t>(std::floor(pos.x() * cellScale));
-                        int64_t by = static_cast<int64_t>(std::floor(pos.y() * cellScale));
-                        int64_t bz = static_cast<int64_t>(std::floor(pos.z() * cellScale));
-
-                        for (int ox = -1; ox <= 1; ox++) {
-                            for (int oy = -1; oy <= 1; oy++) {
-                                for (int oz = -1; oz <= 1; oz++) {
-                                    uint64_t h = ((bx + ox) * 73856093) ^ ((by + oy) * 19349663) ^ ((bz + oz) * 83492791);
-                                    int bucket = h % tableSize;
-                                    
-                                    int curr = head[bucket];
-                                    while (curr != -1) {
-                                        Particle* p = hullPtrs[curr];
-                                        float d2 = (p->currentPos - pos).squaredNorm();
-
-                                        if (d2 < dists[2]) {
-                                            if (d2 < dists[1]) {
-                                                if (d2 < dists[0]) {
-                                                    dists[2] = dists[1];
-                                                    nearest[2] = nearest[1];
-                                                    dists[1] = dists[0];
-                                                    nearest[1] = nearest[0];
-                                                    dists[0] = d2;
-                                                    nearest[0] = p;
-                                                } else {
-                                                    dists[2] = dists[1];
-                                                    nearest[2] = nearest[1];
-                                                    dists[1] = d2;
-                                                    nearest[1] = p;
-                                                }
-                                            } else {
-                                                dists[2] = d2;
-                                                nearest[2] = p;
-                                            }
-                                        }
-                                        curr = next[curr];
-                                    }
-                                }
-                            }
-                        }
-
-                        if (nearest[2] == nullptr) continue;
-
-                        v3 p1 = nearest[0]->currentPos;
-                        v3 p2 = nearest[1]->currentPos;
-                        v3 p3 = nearest[2]->currentPos;
-
-                        v3 v1 = p2 - p1;
-                        v3 v2 = p3 - p1;
-                        v3 normal = v1.cross(v2).normalized();
-
-                        v3 triCenter = (p1 + p2 + p3) * 0.3333f;
-                        if (triCenter.dot(normal) < 0) normal = -normal;
-
-                        v3 toCand = pos - p1;
-                        float dotVal = toCand.dot(normal);
-
-                        if (dotVal < -1e-4f) {
-                            Particle newPt;
-                            newPt.currentPos = pos;
-                            newPt.tectonicPos = pos.cast<Eigen::half>();
-                            newPt.originalPos = pos.cast<Eigen::half>(); 
-                            newPt.surface = false;
-                            
-                            newPt.plateID = nearest[0]->plateID;
-                            newPt.originColor = nearest[0]->originColor;
-
-                            if (std::sqrt(dists[0]) < config.voxelSize * 1.5f) {
-                                tData.surfaceCandidates.push_back(newPt);
-                            } else {
-                                tData.gridCandidates.push_back(newPt);
                         }
                     }
                 }
             }
         }
 
-            #pragma omp critical
-            {
-                config.interpolatedNodes.insert(config.interpolatedNodes.end(), 
-                                                tData.surfaceCandidates.begin(), 
-                                                tData.surfaceCandidates.end());
-                
-                for(const auto& p : tData.surfaceCandidates) {
-                    grid.set(p, p.currentPos, true, p.originColor.cast<float>(), config.voxelSize, true, 1, 0, false, 0.0f, 0.0f, 0.0f);
-                }
-                for(const auto& p : tData.gridCandidates) {
-                    grid.set(p, p.currentPos, true, p.originColor.cast<float>(), config.voxelSize, true, 1, 0, false, 0.0f, 0.0f, 0.0f);
-                }
-            }
-        }
-
-        std::cout << "Volume Fill Complete." << std::endl;
         grid.optimize();
+        std::cout << "Volume Fill Complete. Inserted " << fillCount << " interior nodes directly into the grid." << std::endl;
     }
 
     void simulateImpacts() {