// VoxelMarchingCubesMesher.cpp
// Implémentation du marching cubes density-only.

#include "VoxelMarchingCubesMesher.h"
#include "VoxelHybridSharpMesher.h"
#include "MarchingCubesTables.h"
#include "StandaloneShapeWorld.h"
#include "VoxelForgeCore/Linux77Profiler.h"
#include "VoxelDensityProfile.h"
#if WITH_EDITOR
#include "VoxelRedefinedDensityField.h"
#endif
#include "VoxelStrateManager.h"
#include "VoxelTilePostReach.h"
#include "VoxelPassageGeometry.h"
#include "ProfilingDebugging/CpuProfilerTrace.h"
#include "HAL/IConsoleManager.h"
#include "Misc/CommandLine.h"
#include <atomic>
#include <array>
#include <algorithm>
#include <iterator>
#include <bit>
#include <cstdlib>
#include <immintrin.h>
#include <unordered_map>

namespace
{
    struct FVoxelMesherThreadState
    {
        int32 OctaveBias=0,SampleStep=1,TileCells=0,TileReachBlocks=0,TileReachBlockIndex=INDEX_NONE;
        FIntVector TileOrigin=FIntVector::ZeroValue;
        uint8 TilePostReachFlags=VoxelGenLOD::AllTileReach;
        bool bTilePostReachDebug=false,bTilePostReachBypass=false;
        bool bTileReachDiagnosticTile=false,bTileCoreReachProofEnabled=false;
        bool bTileCoreReachDecisionRecorded=false,bTileBlockReachValid=false;
        float TileReachScale=1.0f;
        std::array<uint8,VoxelGenLOD::TileReachMaxBlockCount> TileBlockReachFlags{};

        static FVoxelMesherThreadState Capture()
        {
            FVoxelMesherThreadState State;
            State.OctaveBias=VoxelGenLOD::OctaveBias;
            State.SampleStep=VoxelGenLOD::SampleStep;
            State.TileOrigin=VoxelGenLOD::TileOriginVoxels;
            State.TileCells=VoxelGenLOD::TileCellsPerAxis;
            State.TilePostReachFlags=VoxelGenLOD::TilePostReachFlags;
            State.bTilePostReachDebug=VoxelGenLOD::bTilePostReachDebug;
            State.bTilePostReachBypass=VoxelGenLOD::bTilePostReachBypass;
            State.bTileReachDiagnosticTile=VoxelGenLOD::bTileReachDiagnosticTile;
            State.bTileCoreReachProofEnabled=VoxelGenLOD::bTileCoreReachProofEnabled;
            State.bTileCoreReachDecisionRecorded=VoxelGenLOD::bTileCoreReachDecisionRecorded;
            State.bTileBlockReachValid=VoxelGenLOD::bTileBlockReachValid;
            State.TileReachScale=VoxelGenLOD::TileReachScale;
            State.TileReachBlocks=VoxelGenLOD::TileReachBlocksPerAxis;
            State.TileReachBlockIndex=VoxelGenLOD::TileReachBlockIndex;
            std::copy(std::begin(VoxelGenLOD::TileBlockReachFlags),
                std::end(VoxelGenLOD::TileBlockReachFlags),State.TileBlockReachFlags.begin());
            return State;
        }

        void Apply() const
        {
            VoxelGenLOD::OctaveBias=OctaveBias;
            VoxelGenLOD::SampleStep=SampleStep;
            VoxelGenLOD::TileOriginVoxels=TileOrigin;
            VoxelGenLOD::TileCellsPerAxis=TileCells;
            VoxelGenLOD::TilePostReachFlags=TilePostReachFlags;
            VoxelGenLOD::bTilePostReachDebug=bTilePostReachDebug;
            VoxelGenLOD::bTilePostReachBypass=bTilePostReachBypass;
            VoxelGenLOD::bTileReachDiagnosticTile=bTileReachDiagnosticTile;
            VoxelGenLOD::bTileCoreReachProofEnabled=bTileCoreReachProofEnabled;
            VoxelGenLOD::bTileCoreReachDecisionRecorded=bTileCoreReachDecisionRecorded;
            VoxelGenLOD::bTileBlockReachValid=bTileBlockReachValid;
            VoxelGenLOD::TileReachScale=TileReachScale;
            VoxelGenLOD::TileReachBlocksPerAxis=TileReachBlocks;
            VoxelGenLOD::TileReachBlockIndex=TileReachBlockIndex;
            std::copy(TileBlockReachFlags.begin(),TileBlockReachFlags.end(),
                std::begin(VoxelGenLOD::TileBlockReachFlags));
        }
    };

    struct FVoxelMesherThreadStateScope
    {
        explicit FVoxelMesherThreadStateScope(const FVoxelMesherThreadState& InTarget)
            :Previous(FVoxelMesherThreadState::Capture())
        { InTarget.Apply(); }
        ~FVoxelMesherThreadStateScope(){Previous.Apply();}
        FVoxelMesherThreadState Previous;
    };

    int32 GVoxelForgeUseOperatorBlock = 1;
    FAutoConsoleVariableRef CVarVoxelForgeUseOperatorBlock(
        TEXT("voxel.UseOperatorBlock"),
        GVoxelForgeUseOperatorBlock,
        TEXT("Use the operator-major density block path for non-shared mesher grids."));
    bool GVoxelForgeBlockSwitchParsed = false;
    bool VF_IsL75TnBlockRegressionEnabled()
    {
        static const bool Enabled=[]
        {
            const char* Value=std::getenv("VOXELFORGE_L75_TN_BLOCK_REGRESSION");
            return Value&&Value[0]=='1'&&Value[1]=='\0';
        }();
        return Enabled;
    }
    int32 GVoxelForgeRawSkirtsForAudit = 0;
    bool GVoxelForgeRawSkirtsSwitchParsed = false;

    void VF_ParseRawSkirtsAuditSwitch()
    {
        if(GVoxelForgeRawSkirtsSwitchParsed)return;
        GVoxelForgeRawSkirtsSwitchParsed=true;
        GVoxelForgeRawSkirtsForAudit=FString(FCommandLine::Get()).Contains(
            TEXT("--audit-raw-skirts"))?1:0;
        if(GVoxelForgeRawSkirtsForAudit)
            std::printf("[SkirtAuditControl] raw skirt curtains enabled for negative-control scan\n");
    }

    int32 GVoxelForgeTilePostReachDebug = 0;
    FAutoConsoleVariableRef CVarVoxelForgeTilePostReachDebug(
        TEXT("voxel.TilePostReachDebug"),
        GVoxelForgeTilePostReachDebug,
        TEXT("Run skipped landing/structural posts in a bypass and compare their exact result; "
             "a nonzero difference fails the debug run."));

    float GVoxelForgeTilePostReachScale = 1.0f;
    FAutoConsoleVariableRef CVarVoxelForgeTilePostReachScale(
        TEXT("voxel.TilePostReachScale"),
        GVoxelForgeTilePostReachScale,
        TEXT("Scale the per-tile post reach for deliberate soundness tests. Default 1.0 is exact; "
             "values below 1 intentionally shrink the guard."));

    int32 GVoxelForgeTileReachCostDiagnostics = 0;
    FAutoConsoleVariableRef CVarVoxelForgeTileReachCostDiagnostics(
        TEXT("voxel.TileReachCostDiagnostics"),
        GVoxelForgeTileReachCostDiagnostics,
        TEXT("Measure core/passage density cost split by tile reach; no field gate."));

    bool GVoxelForgeTilePostReachSwitchesParsed = false;

    void VF_ParseTilePostReachSwitches()
    {
        if (GVoxelForgeTilePostReachSwitchesParsed)
        {
            return;
        }
        GVoxelForgeTilePostReachSwitchesParsed = true;
        int32 DebugValue = GVoxelForgeTilePostReachDebug;
        float ScaleValue = GVoxelForgeTilePostReachScale;
        int32 CostDiagnosticsValue = GVoxelForgeTileReachCostDiagnostics;
        FParse::Value(
            FCommandLine::Get(), TEXT("voxel.TilePostReachDebug="), DebugValue);
        FParse::Value(
            FCommandLine::Get(), TEXT("voxel.TilePostReachScale="), ScaleValue);
        FParse::Value(
            FCommandLine::Get(), TEXT("voxel.TileReachCostDiagnostics="), CostDiagnosticsValue);
        GVoxelForgeTilePostReachDebug = DebugValue != 0 ? 1 : 0;
        if (!VoxelMath::IsFinite(ScaleValue) || ScaleValue <= 0.0f)
        {
            ScaleValue = 1.0f;
        }
        GVoxelForgeTilePostReachScale = FMath::Min(ScaleValue, 1.0f);
        GVoxelForgeTileReachCostDiagnostics = CostDiagnosticsValue != 0 ? 1 : 0;
        VoxelGenLOD::GTilePostReachDebugEnabled.store(
            GVoxelForgeTilePostReachDebug != 0, std::memory_order_relaxed);
        VoxelGenLOD::GTileReachCostDiagnosticsEnabled.store(
            GVoxelForgeTileReachCostDiagnostics != 0, std::memory_order_relaxed);
    }

    void VF_ParseOperatorBlockSwitch()
    {
        if (GVoxelForgeBlockSwitchParsed) { return; }
        GVoxelForgeBlockSwitchParsed = true;
        int32 CommandLineValue = GVoxelForgeUseOperatorBlock;
        if (FParse::Value(FCommandLine::Get(), TEXT("voxel.UseOperatorBlock="), CommandLineValue))
        {
            GVoxelForgeUseOperatorBlock = CommandLineValue;
        }
    }
}

void UVoxelMarchingCubesMesher::InitializeConsoleSwitches()
{
    VF_ParseOperatorBlockSwitch();
    VF_ParseTilePostReachSwitches();
    VF_ParseRawSkirtsAuditSwitch();
}

//=============================================================================
// MAIN ALGORITHM
//=============================================================================
// La grille pré-échantillonnée fournit positions ET gradients.

FVoxelMeshData UVoxelMarchingCubesMesher::GenerateMesh(FIntVector OriginVoxels, int32 Step, int32 InCellsPerAxis,
                                                       TArray<uint8>* OutCaptureGrid,
                                                       int32 BandZMinVox, int32 BandZMaxVox,
                                                       int64* OutDensitySampleCount,
                                                       FVoxelSharedDensityGrid* OutPreEditDensityGrid,
                                                       const FVoxelSharedDensityGrid* InSharedDensityGrid,
                                                       std::vector<float>* OutDebugDensitySamples,
                                                       FVoxelMarchingCubesStageTimings* OutStageTimings,
                                                       FVoxelUniformBlockProofAuditStats* OutBlockProofAudit,
                                                       bool bForceFullDensityCapture)
{
    TRACE_CPUPROFILER_EVENT_SCOPE(VoxelForge_MesherGenerateMesh);
    FStandaloneShapeWorld::FScopedPassageDensityQueryCaller densityQueryCaller(
        EPassageDensityQueryCaller::DensityGrid);
    const uint64 StageMesherStartCycles=OutStageTimings?FPlatformTime::Cycles64():0;
    FVoxelMeshData MeshData;
    if (OutCaptureGrid) { OutCaptureGrid->Reset(); }
    if (OutPreEditDensityGrid) { OutPreEditDensityGrid->Reset(); }
    if (OutDensitySampleCount) { *OutDensitySampleCount = 0; }
    if (OutStageTimings) { *OutStageTimings = FVoxelMarchingCubesStageTimings{}; }
    if (OutBlockProofAudit) { *OutBlockProofAudit = FVoxelUniformBlockProofAuditStats{}; }
    if (OutDebugDensitySamples) { OutDebugDensitySamples->clear(); }
    if (!Generator) return MeshData;
    const bool bEmitUV = !bStandaloneSkipUvAndColors && !bStandaloneSkipUV;
    const bool bEmitColors = !bStandaloneSkipUvAndColors && !bStandaloneSkipColors;
    const bool bUseCoarseFarBiomeGrid = Step >= 64;
    VoxelDensityProfile::FScopedTimer MesherProfileTimer(
        VoxelDensityProfile::EBucket::MesherGenerateMesh);
    VoxelDensityProfile::FScopedTimer MesherSetupTimer(
        VoxelDensityProfile::EBucket::MesherOther);

    // Cell size in voxels. No upper clamp: coarse clipmap levels use bigger steps (the EXTENT
    // grows). Coarse tiles also use FEWER cells (InCellsPerAxis) for cheaper gen.
    Step = FMath::Max(1, Step);
    const float GeometryVoxelSize = VoxelMath::IsFinite(GeometryVoxelSizeCm)
        ? FMath::Max(1.0f, GeometryVoxelSizeCm) : VOXEL_SIZE;

    // T2.b — octave bias for THIS tile: drop LODOctaveDrop octaves per Step doubling from
    // the generator's per-voxel volumetric noise (sub-cell octaves can't shape a coarse
    // isosurface). TGuardValue restores 0 on every exit path, so nothing outside this tile
    // (deco snapping, density-volume fill, the next task on this pooled thread) sees a bias.
    const int32 OctaveBias = (LODOctaveDrop > 0 && Step > 1)
        ? LODOctaveDrop * (int32)FMath::FloorLog2((uint32)Step) : 0;
    TGuardValue<int32> OctaveBiasGuard(VoxelGenLOD::OctaveBias, OctaveBias);
    TGuardValue<int32> SampleStepGuard(VoxelGenLOD::SampleStep, Step);

    // World-cm origin of the tile's min corner (positions are built relative to this).
    const FVector ChunkWorldPos = FVector(OriginVoxels) * GeometryVoxelSize;

    //=========================================================================
    // VERTEX DEDUPLICATION MAP
    //=========================================================================
    // Clé = position quantifiée au 0.01 unité (FIntVector).
    // Valeur = index dans MeshData.Vertices.
    // Les vertices partagés permettent des normales lisses et réduisent le count ~3x.
    // Réutilisé d'une tuile à l'autre (thread_local) : Reset garde les buckets alloués →
    // plus de (ré)allocation de hash-map par tuile (chaque worker a son propre exemplaire).
    static thread_local TMap<FIntVector, int32> VertexMap;
    VertexMap.Reset();
    static thread_local TMap<FIntVector, TMap<uint64, int32>> SplitNormalVertexMap;
    SplitNormalVertexMap.Reset();
    const int32 VertexCapacityHint=FMath::Clamp(InCellsPerAxis,2,CHUNK_SIZE);
    VertexMap.Reserve(VertexCapacityHint*VertexCapacityHint*2);

    //=========================================================================
    // F17 — CLASSE DE SURFACE (sol vs plafond sky-cap), par vertex → par triangle
    //=========================================================================
    // Le discriminant est SÉMANTIQUE, pas géométrique (fable-idea F17) : un vertex orienté vers
    // le bas est un sky-cap seulement s'il est proche de CeilSurf ; proche de TerrainZ c'est un
    // surplomb de terrain (reste "sol" — un vote par tuile mettrait le matériau ciel sous les
    // surplombs des tuiles majoritairement plafond). Un futur toit de grotte (aussi down-
    // facing, mais SOUS TerrainZ) tombera correctement côté "sol/roche" par la même règle.
    // Coût : GetSurfaceHeightAt (pile XY complète) UNIQUEMENT pour les vertex down-facing,
    // mémoïsé par colonne quantifiée au pas de la grille → ≤ (colonnes touchées) appels ;
    // ~zéro sur une tuile de sol pur, borné par lattice² sur une tuile de cap.
    static thread_local TArray<uint8>              VertexClasses;   // 0 = sol, 1 = sky-cap
    static thread_local TMap<FIntVector, FVector2f> SurfColMemo;    // (Xq,Yq,chunkZ) → (TerrainZ, CeilSurf)
    if (!bStandaloneSkipMaterialGroups)
    {
        VertexClasses.Reset();
        SurfColMemo.Reset();
    }

    struct FBiomeGridSample
    {
        int32 Dominant = 0;
        int32 Neighbor = 0;
        float Blend = 0.0f;
    };
    struct FBiomeTileGridSamples
    {
        FBiomeGridSample Corners[4];
    };
    static thread_local TMap<int32, FBiomeTileGridSamples> BiomeTileGridMemo;
    const int32 BiomeGridSpacing = bUseCoarseFarBiomeGrid
        ? FMath::Clamp(InCellsPerAxis, 2, CHUNK_SIZE) * Step : 1;
    if (bUseCoarseFarBiomeGrid && bEmitColors)
    {
        BiomeTileGridMemo.Reset();
    }

    auto GetBiomeTileGridSamples = [&](int32 ChunkZ) -> FBiomeTileGridSamples
    {
        if (const FBiomeTileGridSamples* Existing = BiomeTileGridMemo.Find(ChunkZ))
        {
            return *Existing;
        }
        FBiomeTileGridSamples Samples;
        const float SampleZ = static_cast<float>(ChunkZ * CHUNK_SIZE + CHUNK_SIZE / 2);
        const int32 TileMaxX = OriginVoxels.X + BiomeGridSpacing;
        const int32 TileMaxY = OriginVoxels.Y + BiomeGridSpacing;
        Generator->GetBiomeMaterialAt(static_cast<float>(OriginVoxels.X),
            static_cast<float>(OriginVoxels.Y), SampleZ,
            Samples.Corners[0].Dominant, Samples.Corners[0].Neighbor, Samples.Corners[0].Blend);
        Generator->GetBiomeMaterialAt(static_cast<float>(TileMaxX),
            static_cast<float>(OriginVoxels.Y), SampleZ,
            Samples.Corners[1].Dominant, Samples.Corners[1].Neighbor, Samples.Corners[1].Blend);
        Generator->GetBiomeMaterialAt(static_cast<float>(OriginVoxels.X),
            static_cast<float>(TileMaxY), SampleZ,
            Samples.Corners[2].Dominant, Samples.Corners[2].Neighbor, Samples.Corners[2].Blend);
        Generator->GetBiomeMaterialAt(static_cast<float>(TileMaxX),
            static_cast<float>(TileMaxY), SampleZ,
            Samples.Corners[3].Dominant, Samples.Corners[3].Neighbor, Samples.Corners[3].Blend);
        return BiomeTileGridMemo.Add(ChunkZ, Samples);
    };

    // Normale fournie par l'appelant (gradient lu dans la grille de densité, T1.b).
    auto GetOrCreateVertex = [&](const FVector& WorldPos, const FVector& RawNormal,
                                 uint64 NormalSplitKey = 0) -> int32
    {
        const FIntVector Key(
            FMath::RoundToInt(WorldPos.X * 100.0f),
            FMath::RoundToInt(WorldPos.Y * 100.0f),
            FMath::RoundToInt(WorldPos.Z * 100.0f)
        );

        if (NormalSplitKey != 0)
        {
            if (TMap<uint64, int32>* SplitMap = SplitNormalVertexMap.Find(Key))
            {
                if (int32* Existing = SplitMap->Find(NormalSplitKey)) return *Existing;
            }
        }
        else if (int32* Existing = VertexMap.Find(Key))
        {
            return *Existing;
        }

        const int32 NewIndex = MeshData.Vertices.Num();
        if (NormalSplitKey != 0)
        {
            SplitNormalVertexMap.FindOrAdd(Key).Add(NormalSplitKey, NewIndex);
        }
        else
        {
            VertexMap.Add(Key, NewIndex);
        }

        MeshData.Vertices.Add(WorldPos);

        FVector Normal = RawNormal;
        if (!Normal.Normalize())
        {
            Normal = FVector(0.0f, 0.0f, 1.0f);  // dégénéré (zone plate)
        }
        MeshData.Normals.Add(Normal);

        if (bEmitUV)
        {
            // UVs planaires — le triplanar mapping se fait dans le matériau.
            MeshData.UVs.Add(FVector2f(WorldPos.X / VOXEL_SIZE, WorldPos.Y / VOXEL_SIZE));
        }

        if (bEmitColors)
        {
            //-----------------------------------------------------------------
            // VERTEX COLOUR — masks for biome palette, slope and border blend.
            int32 PalD = 0, PalN = 0; float BlendW = 0.0f;
            if (bUseCoarseFarBiomeGrid)
            {
                const float VoxelX = WorldPos.X / VOXEL_SIZE;
                const float VoxelY = WorldPos.Y / VOXEL_SIZE;
                const float VoxelZ = WorldPos.Z / VOXEL_SIZE;
                const float TX = FMath::Clamp(
                    (VoxelX - static_cast<float>(OriginVoxels.X))
                        / static_cast<float>(BiomeGridSpacing), 0.0f, 1.0f);
                const float TY = FMath::Clamp(
                    (VoxelY - static_cast<float>(OriginVoxels.Y))
                        / static_cast<float>(BiomeGridSpacing), 0.0f, 1.0f);
                const int32 ChunkZ = FMath::FloorToInt(VoxelZ / static_cast<float>(CHUNK_SIZE));
                const FBiomeTileGridSamples TileSamples = GetBiomeTileGridSamples(ChunkZ);
                const FBiomeGridSample& S00 = TileSamples.Corners[0];
                const FBiomeGridSample& S10 = TileSamples.Corners[1];
                const FBiomeGridSample& S01 = TileSamples.Corners[2];
                const FBiomeGridSample& S11 = TileSamples.Corners[3];

                const FBiomeGridSample* Nearest = &S00;
                if (TX >= 0.5f) Nearest = TY >= 0.5f ? &S11 : &S10;
                else if (TY >= 0.5f) Nearest = &S01;
                PalD = Nearest->Dominant;
                PalN = Nearest->Neighbor;
                BlendW = Nearest->Blend;
                if (S00.Dominant == PalD && S00.Neighbor == PalN
                    && S10.Dominant == PalD && S10.Neighbor == PalN
                    && S01.Dominant == PalD && S01.Neighbor == PalN
                    && S11.Dominant == PalD && S11.Neighbor == PalN)
                {
                    const float B0 = FMath::Lerp(S00.Blend, S10.Blend, TX);
                    const float B1 = FMath::Lerp(S01.Blend, S11.Blend, TX);
                    BlendW = FMath::Lerp(B0, B1, TY);
                }
            }
            else
            {
                Generator->GetBiomeMaterialAt(WorldPos.X / VOXEL_SIZE, WorldPos.Y / VOXEL_SIZE,
                                              WorldPos.Z / VOXEL_SIZE, PalD, PalN, BlendW);
            }
            const uint8 R = (uint8)FMath::Clamp(PalD, 0, 255);
            const uint8 A = (uint8)FMath::Clamp(PalN, 0, 255);
            const uint8 G = (uint8)FMath::Clamp(FMath::RoundToInt((1.0f - FMath::Abs((float)Normal.Z)) * 255.0f), 0, 255);
            const uint8 Bc = (uint8)FMath::Clamp(FMath::RoundToInt(BlendW * 255.0f), 0, 255);
            MeshData.Colors.Add(FColor(R, G, Bc, A));
        }

        if (!bStandaloneSkipMaterialGroups)
        {
            // F17 — classify down-facing vertices for material-only ceiling grouping.
            uint8 SurfClass = 0;
            if (Normal.Z < -0.1f)
            {
                const float Zv = WorldPos.Z / GeometryVoxelSize;
                const float StepF = (float)Step;
                const FIntVector MemoKey(
                    FMath::RoundToInt(WorldPos.X / (GeometryVoxelSize * StepF)),
                    FMath::RoundToInt(WorldPos.Y / (GeometryVoxelSize * StepF)),
                    FMath::FloorToInt(Zv / (float)CHUNK_SIZE));
                FVector2f* Col = SurfColMemo.Find(MemoKey);
                if (!Col)
                {
                    float TerrainZ = 0.0f, CeilSurf = 0.0f;
                    bool bSurf=ShapeWorldContext&&ShapeWorldContext->GetNativeSurfaceHeightsAt(
                        static_cast<float>(MemoKey.X*Step),static_cast<float>(MemoKey.Y*Step),
                        MemoKey.Z,TerrainZ,CeilSurf);
                    if(!bSurf)
                    {
                        const float LegacyPerGeometryVoxel = GeometryVoxelSize / VOXEL_SIZE;
                        const float LegacyX = (float)(MemoKey.X * Step) * LegacyPerGeometryVoxel;
                        const float LegacyY = (float)(MemoKey.Y * Step) * LegacyPerGeometryVoxel;
                        const int32 LegacyChunkZ = FMath::FloorToInt(
                            (WorldPos.Z / VOXEL_SIZE) / (float)CHUNK_SIZE);
                        bSurf = Generator->GetSurfaceHeightAt(
                                LegacyX, LegacyY, LegacyChunkZ, TerrainZ, CeilSurf)
                            || Generator->GetSurfaceHeightAt(
                                LegacyX, LegacyY, LegacyChunkZ - 1, TerrainZ, CeilSurf);
                        if (bSurf)
                        {
                            const float LegacyToGeometry = VOXEL_SIZE / GeometryVoxelSize;
                            TerrainZ *= LegacyToGeometry;
                            CeilSurf *= LegacyToGeometry;
                        }
                    }
                    Col = &SurfColMemo.Add(MemoKey, bSurf ? FVector2f(TerrainZ, CeilSurf)
                                                          : FVector2f(FLT_MAX, -FLT_MAX));
                }
                if (Col->X != FLT_MAX && FMath::Abs(Zv - Col->Y) < FMath::Abs(Zv - Col->X)
                    && Zv <= Col->Y + 2.0f * StepF)
                {
                    SurfClass = 1;
                }
            }
            VertexClasses.Add(SurfClass);
        }

        return NewIndex;
    };

    //=========================================================================
    // CORNER + EDGE TABLES
    //=========================================================================
    //       4-------5
    //      /|      /|
    //     / |     / |
    //    7-------6  |
    //    |  0----|--1
    //    | /     | /
    //    |/      |/
    //    3-------2
    static const FIntVector CornerOffsets[8] = {
        {0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0},
        {0, 0, 1}, {1, 0, 1}, {1, 1, 1}, {0, 1, 1},
    };

    static const int32 EdgeCorners[12][2] = {
        {0, 1}, {1, 2}, {2, 3}, {3, 0},  // Bas (0-3)
        {4, 5}, {5, 6}, {6, 7}, {7, 4},  // Haut (4-7)
        {0, 4}, {1, 5}, {2, 6}, {3, 7},  // Verticales (8-11)
    };

    //=========================================================================
    // PRÉ-CALCUL DE LA GRILLE DE DENSITÉ
    //=========================================================================
    // Les cellules adjacentes partagent leurs coins : en échantillonnant par
    // cellule (8 coins) on appelle GetDensityAt ~8× de trop pour chaque point.
    // On échantillonne donc chaque point de grille UNE SEULE FOIS dans un tableau
    // plat, puis le balayage des cellules y lit ses coins. GetDensityAt est une
    // fonction pure de la coordonnée monde, donc le maillage est identique au bit
    // près — c'est juste ~7× moins d'appels au LOD0 (33³ au lieu de 32³×8).
    //
    // CellsPerAxis = nombre de cellules par axe ; GridDim = points de grille (+1).
    // Le point de grille (gx,gy,gz) correspond au voxel monde (gx,gy,gz)*Step.
    // The grid fill order is selected below by tile extent to preserve room-cache locality.
    // On échantillonne avec un anneau de marge de 1 point de chaque côté (indices -1..GridDim)
    // pour pouvoir calculer les normales par GRADIENT DE GRILLE (T1.b) — différences centrales
    // sur la grille au lieu de 6 appels densité frais par vertex. La marge utilise les mêmes
    // échantillons monde purs qu'un chunk voisin, donc les normales restent continues aux bords
    // de chunk (pas de couture de shading). La géométrie est inchangée au bit près (mêmes
    // positions d'arête) ; seules les normales changent.
    const int32 CellsPerAxis = FMath::Clamp(InCellsPerAxis, 2, CHUNK_SIZE);  // coarse tiles use fewer
    const int32 GridDim      = CellsPerAxis + 1;
    const int32 MDim         = GridDim + 2;            // +1 marge de chaque côté
    TGuardValue<FIntVector> TileOriginGuard(
        VoxelGenLOD::TileOriginVoxels, OriginVoxels);
    TGuardValue<int32> TileCellsGuard(
        VoxelGenLOD::TileCellsPerAxis, CellsPerAxis);

    // Every density sample in this tile sees the same post-reach decision.  The one-point halo
    // matches the actual grid fill below (gx/gy/gz = -1..GridDim), so a false bit proves that no
    // landing or passage structural writer can reach any sample used by meshing or gradients.
    TGuardValue<bool> TilePostReachDebugGuard(
        VoxelGenLOD::bTilePostReachDebug,
        GVoxelForgeTilePostReachDebug != 0);
    TGuardValue<bool> TilePostReachBypassGuard(
        VoxelGenLOD::bTilePostReachBypass, false);
    TGuardValue<bool> TileReachDiagnosticGuard(
        VoxelGenLOD::bTileReachDiagnosticTile, true);
    TGuardValue<bool> TileCoreReachProofGuard(
        VoxelGenLOD::bTileCoreReachProofEnabled, false);
    TGuardValue<bool> TileCoreReachDecisionGuard(
        VoxelGenLOD::bTileCoreReachDecisionRecorded, false);
    TGuardValue<float> TileReachScaleGuard(
        VoxelGenLOD::TileReachScale, GVoxelForgeTilePostReachScale);
    const float ReachScale = GVoxelForgeTilePostReachScale;
    const int32 ReachBlocksPerAxis = FMath::DivideAndRoundUp(
        CellsPerAxis, VoxelGenLOD::TileReachBlockCells);
    TGuardValue<int32> TileReachBlocksGuard(
        VoxelGenLOD::TileReachBlocksPerAxis, ReachBlocksPerAxis);
    TGuardValue<int32> TileReachBlockIndexGuard(
        VoxelGenLOD::TileReachBlockIndex, INDEX_NONE);
    TGuardValue<bool> TileBlockReachValidGuard(
        VoxelGenLOD::bTileBlockReachValid, false);
    uint8 TilePostReachFlags = VoxelGenLOD::AllTilePostReach;
    if (Generator->StrateManager != nullptr)
    {
        const int64 Extent = static_cast<int64>(CellsPerAxis + 1) * Step;
        const FBox LatticeWithHalo(
            FVector(
                static_cast<float>(static_cast<int64>(OriginVoxels.X) - Step),
                static_cast<float>(static_cast<int64>(OriginVoxels.Y) - Step),
                static_cast<float>(static_cast<int64>(OriginVoxels.Z) - Step)),
            FVector(
                static_cast<float>(static_cast<int64>(OriginVoxels.X) + Extent),
                static_cast<float>(static_cast<int64>(OriginVoxels.Y) + Extent),
                static_cast<float>(static_cast<int64>(OriginVoxels.Z) + Extent)));
        // Core stays conservative until the tile-window tunnel cache is built.  Passage carving
        // is manager-owned and can be decided immediately from the same lattice+halo box.
        TilePostReachFlags = VoxelGenLOD::TunnelCoreReachable;
        if (Generator->StrateManager->AnyOriginLandingNearLattice(
                LatticeWithHalo, OriginVoxels, Step, ReachScale))
        {
            TilePostReachFlags |= VoxelGenLOD::OriginLandingReachable;
        }
        if (Generator->StrateManager->AnyPassageLandingNearLattice(
                LatticeWithHalo, OriginVoxels, Step, ReachScale))
        {
            TilePostReachFlags |= VoxelGenLOD::PassageLandingReachable;
        }
        if (Generator->StrateManager->AnyPassageCarvingNearLattice(
                LatticeWithHalo, OriginVoxels, Step, ReachScale))
        {
            TilePostReachFlags |= VoxelGenLOD::PassageCarvingReachable;
        }
        if (Generator->StrateManager->AnyPassageStructuralPostNearLattice(
                LatticeWithHalo, OriginVoxels, Step, ReachScale))
        {
            TilePostReachFlags |= VoxelGenLOD::PassageStructuralReachable;
        }
    }
    TGuardValue<uint8> TilePostReachFlagsGuard(
        VoxelGenLOD::TilePostReachFlags, TilePostReachFlags);
    for (int32 BlockZ = 0; BlockZ < ReachBlocksPerAxis; ++BlockZ)
    {
        for (int32 BlockY = 0; BlockY < ReachBlocksPerAxis; ++BlockY)
        {
            for (int32 BlockX = 0; BlockX < ReachBlocksPerAxis; ++BlockX)
            {
                const int32 BlockIndex =
                    ((BlockZ * ReachBlocksPerAxis) + BlockY) * ReachBlocksPerAxis + BlockX;
                uint8 BlockFlags = TilePostReachFlags;
                if (Generator->StrateManager != nullptr
                    && (TilePostReachFlags & VoxelGenLOD::PassageCarvingReachable) != 0
                    && !Generator->StrateManager->AnyPassageCarvingNearLattice(
                        VoxelGenLOD::MakeTileReachBlockBox(
                            OriginVoxels, Step, CellsPerAxis,
                            BlockX, BlockY, BlockZ),
                        OriginVoxels, Step, ReachScale))
                {
                    BlockFlags &= ~VoxelGenLOD::PassageCarvingReachable;
                }
                VoxelGenLOD::TileBlockReachFlags[BlockIndex] = BlockFlags;
            }
        }
    }
    VoxelGenLOD::bTileBlockReachValid = true;
    VoxelGenLOD::RecordTileReachMask(TilePostReachFlags);
    // COUPE DE CONTENU PAR STRATE — restreint le maillage (et l'échantillonnage) aux cellules
    // dont l'intervalle Z chevauche la bande [BandZMinVox, BandZMaxVox] (voxels inclusifs).
    // Les tuiles à capture ne sont jamais bandées (niveau 0 — garde-fou ci-dessous).
    if (OutCaptureGrid || OutPreEditDensityGrid)
    {
        BandZMinVox = INT32_MIN;
        BandZMaxVox = INT32_MAX;
    }
    int32 CzLo = 0, CzHi = CellsPerAxis - 1;
    if (BandZMinVox > INT32_MIN || BandZMaxVox < INT32_MAX)
    {
        auto FloorDivI = [](int64 A, int64 B) -> int32
        {
            const int64 Q = A / B;
            return (int32)(Q - (((A % B) != 0 && ((A < 0) != (B < 0))) ? 1 : 0));
        };
        // Cellule c couvre [O.Z + c*Step, O.Z + (c+1)*Step] : première/dernière cellule contenant
        // la borne (un contact purement tangent est exclu — sa traversée serait hors bande).
        CzLo = FMath::Max(CzLo, FloorDivI((int64)BandZMinVox - OriginVoxels.Z, Step));
        CzHi = FMath::Min(CzHi, FloorDivI((int64)BandZMaxVox - OriginVoxels.Z, Step));
        if (CzHi < CzLo)
        {
            if(OutStageTimings)
            {
                const uint64 Candidate=static_cast<uint64>(CellsPerAxis)*CellsPerAxis*CellsPerAxis;
                OutStageTimings->CandidateCells=Candidate;
                OutStageTimings->ZBandRejectedCells=Candidate;
            }
            return MeshData;
        }          // tuile entièrement hors bande → vide
    }
    const int32 GzLo = CzLo - 1;                       // coins CzLo..CzHi+1, gradients ±1
    const int32 GzHi = CzHi + 2;                       // (== -1..GridDim sans bande)

    MesherSetupTimer.End();

    // Optional 8x8x8-cell block proof. Mixed/Unknown blocks always fall back to exact density
    // sampling. In the standalone shape route, certified solid/air blocks also omit density-grid
    // fills outside active blocks and their one-sample normal halo.
    TArray<uint8> BlockSkip;
    int32 BlocksPerAxis = 0;
    // Step 2 is enabled for standalone shape proofs (or explicit proof audits). Steps 16/32 keep
    // the historical gate off because their smaller grids do not amortize per-block queries.
    const bool bNativeSurfaceTileProof = Step == 2 && ShapeWorldContext != nullptr
        && ShapeWorldContext->IsNativeSurfaceBox(OriginVoxels,Step,CellsPerAxis);
    const bool bStandaloneShapeProof = ShapeWorldContext != nullptr;
    const bool bUseProofBlocks = bUseBlockEarlyOut && !OutCaptureGrid
        && (Step == 1 || Step == 8 || (Step == 2
            && (bNativeSurfaceTileProof || OutBlockProofAudit != nullptr
                || bStandaloneShapeProof)));
    bool bCanSkipProofBlockGrid = Step <= 2 && bStandaloneShapeProof;
    const int32 ProofBlockCells=std::clamp(bStandaloneProofBlockCells,1,CellsPerAxis);
    if (bUseProofBlocks)
    {
        const int32 BlockCells = ProofBlockCells;
        BlocksPerAxis = (CellsPerAxis + BlockCells - 1) / BlockCells;
        BlockSkip.Init(0, BlocksPerAxis * BlocksPerAxis * BlocksPerAxis);
        bool bHaveUniformBlock = false;
        EVoxelTileClass FirstUniformVerdict = EVoxelTileClass::Mixed;
        bool bAllRelevantBlocksUniform = true;
        bool bAllRelevantBlocksSame = true;
        bool bAllRelevantBlocksProofSafe = true;
        TArray<int32> MixedBlockIndices;

        for (int32 BlockZ = 0; BlockZ < BlocksPerAxis; ++BlockZ)
        {
            const int32 BlockMinZ = BlockZ * BlockCells;
            const int32 BlockMaxZ = FMath::Min(CellsPerAxis - 1,
                                               BlockMinZ + BlockCells - 1);
            if (BlockMaxZ < CzLo || BlockMinZ > CzHi)
            {
                continue;
            }
            for (int32 BlockY = 0; BlockY < BlocksPerAxis; ++BlockY)
            {
                for (int32 BlockX = 0; BlockX < BlocksPerAxis; ++BlockX)
                {
                    const FIntVector BlockOrigin = OriginVoxels + FIntVector(
                        BlockX * BlockCells * Step,
                        BlockY * BlockCells * Step,
                        BlockMinZ * Step);
                    if (VoxelDensityProfile::AreCountersEnabled())
                    {
                        VoxelDensityProfile::AddCounter(
                            VoxelDensityProfile::ECounter::MesherBlockTests);
                    }
                    const bool bAttributeClassifier =
                        VoxelDensityProfile::GetMode()
                            == VoxelDensityProfile::EMode::Attribution;
                    const bool bSampleClassifier = bAttributeClassifier
                        && VoxelDensityProfile::ShouldSample(
                            VoxelDensityProfile::EBucket::ClassifierTotal);
                    const uint64 ClassifierStartCycles = bSampleClassifier
                        ? FPlatformTime::Cycles64() : 0;
                    FVoxelTileClassificationStats ClassifierStats;
                    if (bAttributeClassifier)
                    {
                        VoxelDensityProfile::RecordClassifierCall(bSampleClassifier);
                    }
                    EVoxelTileClass Verdict;
                    bool bBlockUnknown = false;
                    bool bBlockNative = false;
                    bool bBlockProofSafe = false;
                    if(ShapeWorldContext&&Step<=2)
                    {
                        const auto ShapeVerdict=ShapeWorldContext->ClassifyBoxDetailed(
                            BlockOrigin,Step,BlockCells,Step==2);
                        Verdict=!ShapeVerdict.bUnknown
                            &&ShapeVerdict.Verdict==EVoxelTileClass::AllSolid
                            ?EVoxelTileClass::AllSolid
                            :!ShapeVerdict.bUnknown
                                &&ShapeVerdict.Verdict==EVoxelTileClass::AllAir
                                ?EVoxelTileClass::AllAir:EVoxelTileClass::Mixed;
                        bBlockUnknown=ShapeVerdict.bUnknown;
                        bBlockNative=ShapeVerdict.bNativeSurface;
                        bBlockProofSafe=ShapeVerdict.bSkipSafe;
                        if(VF_IsL75TnBlockRegressionEnabled())
                        {
                            constexpr int32 ProbeX=-80,ProbeY=-20,ProbeZ=-3218;
                            const auto IsCapturedGridPoint=[&](int32 Value,int32 Start)
                            {
                                const int32 Delta=Value-Start;
                                return Delta%Step==0&&Delta/Step>=-1
                                    &&Delta/Step<=BlockCells+1;
                            };
                            if(IsCapturedGridPoint(ProbeX,BlockOrigin.X)
                                &&IsCapturedGridPoint(ProbeY,BlockOrigin.Y)
                                &&IsCapturedGridPoint(ProbeZ,BlockOrigin.Z))
                            {
                                const float Density=ShapeWorldContext->DensityAtLattice(
                                    ProbeX,ProbeY,ProbeZ);
                                const std::string Dispatch=ShapeWorldContext->DensityDispatchPathAt(
                                    static_cast<float>(ProbeX),static_cast<float>(ProbeY),
                                    static_cast<float>(ProbeZ));
                                const EVoxelTileClass GeneratorVerdict=Generator->ClassifyTile(
                                    BlockOrigin,Step,BlockCells);
                                const char* Label=Verdict==EVoxelTileClass::AllSolid?"solid":
                                    Verdict==EVoxelTileClass::AllAir?"air":"mixed";
                                const char* GeneratorLabel=GeneratorVerdict==EVoxelTileClass::AllSolid
                                    ?"solid":GeneratorVerdict==EVoxelTileClass::AllAir?"air":"mixed";
                                std::printf("[L75TNBlockRegression] tile_origin=(%d,%d,%d) step=%d cells=%d block=(%d,%d,%d) verdict=%s proof_safe=%d unknown=%d density=%.9g dispatch=%s generator_verdict=%s\n",
                                    OriginVoxels.X,OriginVoxels.Y,OriginVoxels.Z,Step,
                                    CellsPerAxis,BlockX,BlockY,BlockZ,Label,
                                    bBlockProofSafe?1:0,bBlockUnknown?1:0,Density,
                                    Dispatch.c_str(),GeneratorLabel);
                            }
                        }
                        if((Verdict==EVoxelTileClass::AllSolid
                                ||Verdict==EVoxelTileClass::AllAir)
                            &&!bBlockProofSafe)
                        {
                            // An uncertified uniform candidate makes the parent block grid
                            // unusable as a sampling certificate. Keep the committed fallback:
                            // fill the full density grid for this tile.
                            bCanSkipProofBlockGrid=false;
                            bAllRelevantBlocksProofSafe=false;
                        }
                    }
                    else
                    {
                        if(Step==2)bCanSkipProofBlockGrid=false;
                        Verdict=Generator->ClassifyTile(BlockOrigin,Step,BlockCells,
                            bSampleClassifier ? &ClassifierStats : nullptr);
                    }
                    if (bSampleClassifier)
                    {
                        VoxelDensityProfile::AddSampledMeasurement(
                            VoxelDensityProfile::EBucket::ClassifierTotal,
                            FPlatformTime::Cycles64() - ClassifierStartCycles);
                        VoxelDensityProfile::FSnapshot::FClassifierStats ProfileStats;
                        ProfileStats.RefineNodes = ClassifierStats.RefineNodes;
                        ProfileStats.StackBoxCalls = ClassifierStats.StackBoxCalls;
                        ProfileStats.WholeMixedNodes = ClassifierStats.WholeMixedNodes;
                        ProfileStats.WholeSolidNodes = ClassifierStats.WholeSolidNodes;
                        ProfileStats.WholeAirNodes = ClassifierStats.WholeAirNodes;
                        ProfileStats.NeedsFinalFieldNodes = ClassifierStats.NeedsFinalFieldNodes;
                        ProfileStats.SplitNodes = ClassifierStats.SplitNodes;
                        ProfileStats.MaxRefinementDepth = ClassifierStats.MaxRefinementDepth;
                        ProfileStats.ExactCoreSamples = ClassifierStats.ExactCoreSamples;
                        ProfileStats.ExactFinalSamples = ClassifierStats.ExactFinalSamples;
                        ProfileStats.ExactCoreCacheHits = ClassifierStats.ExactCoreCacheHits;
                        ProfileStats.ExactFinalCacheHits = ClassifierStats.ExactFinalCacheHits;
                        ProfileStats.ExactCoreLeaves = ClassifierStats.ExactCoreLeaves;
                        ProfileStats.ExactFinalLeaves = ClassifierStats.ExactFinalLeaves;
                        ProfileStats.StackBoxCycles = ClassifierStats.StackBoxCycles;
                        ProfileStats.ExactCoreCycles = ClassifierStats.ExactCoreCycles;
                        ProfileStats.ExactFinalCycles = ClassifierStats.ExactFinalCycles;
                        ProfileStats.RoomTailQueries = ClassifierStats.RoomTailQueries;
                        ProfileStats.RoomTailEvaluated = ClassifierStats.RoomTailEvaluated;
                        ProfileStats.RoomTailCycles = ClassifierStats.RoomTailCycles;
                        ProfileStats.RoomPropagateCycles = ClassifierStats.RoomPropagateCycles;
                        ProfileStats.RoomExactPrimitiveCycles = ClassifierStats.RoomExactPrimitiveCycles;
                        ProfileStats.RoomCacheWindowCycles = ClassifierStats.RoomCacheWindowCycles;
                        ProfileStats.RoomNumRooms = FMath::Max(ClassifierStats.RoomNumRooms, 0);
                        ProfileStats.RoomNumTunnels = FMath::Max(ClassifierStats.RoomNumTunnels, 0);
                        ProfileStats.RoomNumRoomFloorJoins = FMath::Max(
                            ClassifierStats.RoomNumRoomFloorJoins, 0);
                        ProfileStats.RoomNumPits = FMath::Max(ClassifierStats.RoomNumPits, 0);
                        ProfileStats.RoomNumChimneys = FMath::Max(ClassifierStats.RoomNumChimneys, 0);
                        VoxelDensityProfile::AddClassifierStats(ProfileStats);

                        if (ClassifierStats.StackBoxCalls > 0)
                        {
                            VoxelDensityProfile::AddMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierIntervalProof,
                                0, ClassifierStats.StackBoxCalls);
                            VoxelDensityProfile::AddSampledMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierIntervalProof,
                                ClassifierStats.StackBoxCycles);
                        }
                        if (ClassifierStats.ExactCoreLeaves > 0)
                        {
                            VoxelDensityProfile::AddMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierExactCore,
                                0, ClassifierStats.ExactCoreLeaves);
                            VoxelDensityProfile::AddSampledMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierExactCore,
                                ClassifierStats.ExactCoreCycles);
                        }
                        if (ClassifierStats.ExactFinalLeaves > 0)
                        {
                            VoxelDensityProfile::AddMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierExactFinal,
                                0, ClassifierStats.ExactFinalLeaves);
                            VoxelDensityProfile::AddSampledMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierExactFinal,
                                ClassifierStats.ExactFinalCycles);
                        }
                        if (ClassifierStats.RoomTailQueries > 0)
                        {
                            VoxelDensityProfile::AddMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierRoomTail,
                                0, ClassifierStats.RoomTailQueries);
                            VoxelDensityProfile::AddSampledMeasurement(
                                VoxelDensityProfile::EBucket::ClassifierRoomTail,
                                ClassifierStats.RoomTailCycles);
                        }
                    }
                    if (Verdict == EVoxelTileClass::AllSolid
                        || Verdict == EVoxelTileClass::AllAir)
                    {
                        if(OutBlockProofAudit&&Step<=2&&bBlockProofSafe)
                        {
                            ++OutBlockProofAudit->CheckedBlocks;
                            const int32 ProofMinX=BlockX*BlockCells-1;
                            const int32 ProofMinY=BlockY*BlockCells-1;
                            const int32 ProofMinZ=BlockZ*BlockCells-1;
                            const int32 ProofMaxX=BlockX*BlockCells+BlockCells+1;
                            const int32 ProofMaxY=BlockY*BlockCells+BlockCells+1;
                            const int32 ProofMaxZ=BlockZ*BlockCells+BlockCells+1;
                            bool bUniform=true;
                            for(int32 Z=ProofMinZ;Z<=ProofMaxZ&&bUniform;++Z)
                            for(int32 Y=ProofMinY;Y<=ProofMaxY&&bUniform;++Y)
                            for(int32 X=ProofMinX;X<=ProofMaxX;++X)
                            {
                                ++OutBlockProofAudit->CheckedSamples;
                                const float Density=ShapeWorldContext
                                    ?ShapeWorldContext->DensityAtLattice(
                                        OriginVoxels.X+X*Step,
                                        OriginVoxels.Y+Y*Step,
                                        OriginVoxels.Z+Z*Step)
                                    :Generator->GetDensityAtLattice(
                                        OriginVoxels.X+X*Step,
                                        OriginVoxels.Y+Y*Step,
                                        OriginVoxels.Z+Z*Step);
                                if(!VoxelMath::IsFinite(Density)||Density==0.0f
                                    ||(Verdict==EVoxelTileClass::AllSolid
                                        ?Density>=0.0f:Density<=0.0f))
                                {
                                    std::printf("[UniformBlockProofViolation] tile_origin=(%d,%d,%d) step=%d cells=%d block=(%d,%d,%d) sample=(%d,%d,%d) expected=%s density=%.9g native=%d unknown=%d\n",
                                        OriginVoxels.X,OriginVoxels.Y,OriginVoxels.Z,Step,
                                        CellsPerAxis,BlockX,BlockY,BlockZ,
                                        OriginVoxels.X+X*Step,OriginVoxels.Y+Y*Step,
                                        OriginVoxels.Z+Z*Step,
                                        Verdict==EVoxelTileClass::AllSolid?"solid":"air",
                                        Density,bBlockNative?1:0,bBlockUnknown?1:0);
                                    bUniform=false;
                                    break;
                                }
                            }
                            if(!bUniform)
                            {
                                ++OutBlockProofAudit->Violations;
                                Verdict=EVoxelTileClass::Mixed;
                            }
                        }
                    }
                    if (Verdict == EVoxelTileClass::Mixed)
                    {
                        if(OutBlockProofAudit&&Step<=2)
                        {
                            ++OutBlockProofAudit->FallbackBlocks;
                            if(bBlockUnknown)++OutBlockProofAudit->UnknownBlocks;
                        }
                        MixedBlockIndices.Add(
                            ((BlockZ * BlocksPerAxis) + BlockY) * BlocksPerAxis + BlockX);
                    }
                    else if(OutBlockProofAudit&&Step<=2&&!bBlockProofSafe)
                    {
                        // A midpoint or other non-interval verdict is only a candidate;
                        // it is not counted as a certified hit in the statistical sweep.
                        ++OutBlockProofAudit->FallbackBlocks;
                    }
                    if (Verdict == EVoxelTileClass::AllSolid
                        || Verdict == EVoxelTileClass::AllAir)
                    {
                        // In the ordinary path, candidates are checked against the filled grid.
                        // The Step-2 supremum interval path is itself the proof used to omit fills.
                        BlockSkip[((BlockZ * BlocksPerAxis) + BlockY) * BlocksPerAxis + BlockX] =
                            Verdict == EVoxelTileClass::AllSolid ? 1 : 2;
                        if (!bHaveUniformBlock)
                        {
                            FirstUniformVerdict = Verdict;
                            bHaveUniformBlock = true;
                        }
                        else if (FirstUniformVerdict != Verdict)
                        {
                            bAllRelevantBlocksSame = false;
                        }
                    }
                    else
                    {
                        bAllRelevantBlocksUniform = false;
                    }
                }
            }
        }

        // A conservative whole-tile proof can leave a small block unresolved even when every
        // other block has the same verdict. Resolve that tail on the exact MC lattice. This uses
        // the final GetDensityAt field, checks every cell corner, treats zero/non-finite values
        // as unknown, and never changes the grid resolution or skips an uncertified block.
        if (bHaveUniformBlock && bAllRelevantBlocksSame && !bAllRelevantBlocksUniform
            && MixedBlockIndices.Num() <= 4)
        {
            int32 RemainingMixedBlocks = MixedBlockIndices.Num();
            for (const int32 MixedBlockIndex : MixedBlockIndices)
            {
                const int32 BlockX = MixedBlockIndex % BlocksPerAxis;
                const int32 BlockY = (MixedBlockIndex / BlocksPerAxis) % BlocksPerAxis;
                const int32 BlockZ = MixedBlockIndex / (BlocksPerAxis * BlocksPerAxis);
                const FIntVector BlockOrigin = OriginVoxels + FIntVector(
                    BlockX * BlockCells * Step,
                    BlockY * BlockCells * Step,
                    BlockZ * BlockCells * Step);

                bool bAllSolid = true;
                bool bAllAir = true;
                for (int32 Z = 0; Z <= BlockCells && (bAllSolid || bAllAir); ++Z)
                {
                    for (int32 Y = 0; Y <= BlockCells && (bAllSolid || bAllAir); ++Y)
                    {
                        for (int32 X = 0; X <= BlockCells; ++X)
                        {
                            const float Density = ShapeWorldContext
                                ?ShapeWorldContext->DensityAt(
                                    static_cast<float>(BlockOrigin.X+X*Step),
                                    static_cast<float>(BlockOrigin.Y+Y*Step),
                                    static_cast<float>(BlockOrigin.Z+Z*Step))
                                :Generator->GetDensityAt(BlockOrigin.X+X*Step,
                                    BlockOrigin.Y+Y*Step,BlockOrigin.Z+Z*Step);
                            if (!VoxelMath::IsFinite(Density) || Density == 0.0f)
                            {
                                bAllSolid = false;
                                bAllAir = false;
                                break;
                            }
                            bAllSolid &= Density < 0.0f;
                            bAllAir &= Density > 0.0f;
                            if (!bAllSolid && !bAllAir)
                            {
                                break;
                            }
                        }
                    }
                }

                const EVoxelTileClass ExactVerdict = bAllSolid
                    ? EVoxelTileClass::AllSolid
                    : (bAllAir ? EVoxelTileClass::AllAir : EVoxelTileClass::Mixed);
                if (ExactVerdict == EVoxelTileClass::AllSolid
                    || ExactVerdict == EVoxelTileClass::AllAir)
                {
                    BlockSkip[MixedBlockIndex] =
                        ExactVerdict == EVoxelTileClass::AllSolid ? 1 : 2;
                    --RemainingMixedBlocks;
                    if (ExactVerdict != FirstUniformVerdict)
                    {
                        bAllRelevantBlocksSame = false;
                    }
                }
            }
            bAllRelevantBlocksUniform = RemainingMixedBlocks == 0;
        }

        // Every relevant block has independently proved the same sign.  This is a complete
        // certificate for the tile's MC lattice: the block boxes share their boundary vertices,
        // so equal signs on both sides cannot hide a surface between blocks.  The whole-tile
        // ClassifyTile path already relies on this same non-Mixed proof contract. Return before
        // allocating the density grid; no samples are skipped for a tile that can contain
        // geometry, and the final geometry remains byte-identical because this path emits none.
        if (bHaveUniformBlock && bAllRelevantBlocksSame && bAllRelevantBlocksUniform
            && bAllRelevantBlocksProofSafe
            && !(Step == 2 && OutDebugDensitySamples != nullptr)
            && !(bForceFullDensityCapture && OutDebugDensitySamples != nullptr))
        {
            return MeshData;
        }

        // A mixed or differently signed block set still needs the exact mesher grid.  Classifier
        // results remain candidates only in that case; the grid validation below keeps the cell
        // skip independent of worker-local cache state and preserves geometry byte-for-byte.
    }

    // Réutilise le tampon entre tuiles (thread_local) : SetNumUninitialized garde la
    // capacité, donc plus de malloc/free de ~170 Ko (35³ floats) par tuile.
    static thread_local TArray<float> DensityGrid;
    DensityGrid.SetNumUninitialized(MDim * MDim * MDim);
    float* const DensityGridData=DensityGrid.GetData();
    const FVoxelSharedDensityGrid* ReuseGrid = InSharedDensityGrid != nullptr
        ? InSharedDensityGrid : SharedDensityGrid;
    const int32 SharedDeltaX = ReuseGrid
        ? (OriginVoxels.X - ReuseGrid->OriginVoxels.X) / Step : 0;
    const int32 SharedDeltaY = ReuseGrid
        ? (OriginVoxels.Y - ReuseGrid->OriginVoxels.Y) / Step : 0;
    const int32 SharedDeltaZ = ReuseGrid
        ? (OriginVoxels.Z - ReuseGrid->OriginVoxels.Z) / Step : 0;
    const bool bUseSharedDensityGrid = ReuseGrid
        && ReuseGrid->Step == Step
        && ReuseGrid->Dim > 0
        && ReuseGrid->Samples.Num() == ReuseGrid->Dim * ReuseGrid->Dim * ReuseGrid->Dim
        && (OriginVoxels.X - ReuseGrid->OriginVoxels.X) % Step == 0
        && (OriginVoxels.Y - ReuseGrid->OriginVoxels.Y) % Step == 0
        && (OriginVoxels.Z - ReuseGrid->OriginVoxels.Z) % Step == 0
        && SharedDeltaX - 1 >= 0
        && SharedDeltaY - 1 >= 0
        && SharedDeltaZ - 1 >= 0
        && SharedDeltaX + GridDim < ReuseGrid->Dim
        && SharedDeltaY + GridDim < ReuseGrid->Dim
        && SharedDeltaZ + GridDim < ReuseGrid->Dim;
    const bool bSkipProofBlockGrid = bUseProofBlocks && bCanSkipProofBlockGrid
        && OutBlockProofAudit == nullptr && OutPreEditDensityGrid == nullptr
        && !bUseSharedDensityGrid;
    static thread_local TArray<uint8> RequiredGridSamples;
    if(bSkipProofBlockGrid)
    {
        RequiredGridSamples.SetNumUninitialized(MDim*MDim*MDim);
        FMemory::Memzero(RequiredGridSamples.GetData(),
            static_cast<SIZE_T>(RequiredGridSamples.Num())*sizeof(uint8));
        const int32 BlockCells=ProofBlockCells;
        for(int32 BlockZ=0;BlockZ<BlocksPerAxis;++BlockZ)
        {
            const int32 BlockMinZ=BlockZ*BlockCells;
            const int32 BlockMaxZ=FMath::Min(CellsPerAxis-1,BlockMinZ+BlockCells-1);
            if(BlockMaxZ<CzLo||BlockMinZ>CzHi)continue;
            for(int32 BlockY=0;BlockY<BlocksPerAxis;++BlockY)
            for(int32 BlockX=0;BlockX<BlocksPerAxis;++BlockX)
            {
                const int32 BlockIndex=((BlockZ*BlocksPerAxis)+BlockY)*BlocksPerAxis+BlockX;
                if(BlockSkip[BlockIndex]!=0)continue;
                const int32 CellMinX=BlockX*BlockCells;
                const int32 CellMinY=BlockY*BlockCells;
                const int32 CellMaxX=FMath::Min(CellsPerAxis-1,CellMinX+BlockCells-1);
                const int32 CellMaxY=FMath::Min(CellsPerAxis-1,CellMinY+BlockCells-1);
                const int32 SampleMinX=FMath::Max(-1,CellMinX-1);
                const int32 SampleMaxX=FMath::Min(GridDim,CellMaxX+2);
                const int32 SampleMinY=FMath::Max(-1,CellMinY-1);
                const int32 SampleMaxY=FMath::Min(GridDim,CellMaxY+2);
                const int32 SampleMinZ=FMath::Max(GzLo,FMath::Max(CzLo,BlockMinZ)-1);
                const int32 SampleMaxZ=FMath::Min(GzHi,FMath::Min(CzHi,BlockMaxZ)+2);
                for(int32 Z=SampleMinZ;Z<=SampleMaxZ;++Z)
                for(int32 Y=SampleMinY;Y<=SampleMaxY;++Y)
                for(int32 X=SampleMinX;X<=SampleMaxX;++X)
                    RequiredGridSamples[((Z+1)*MDim+(Y+1))*MDim+(X+1)]=1;
            }
        }
    }
    if(OutStageTimings)
    {
        const uint64 PlaneCells=static_cast<uint64>(CellsPerAxis)*CellsPerAxis;
        const int32 ActiveZ=FMath::Max(0,CzHi-CzLo+1);
        OutStageTimings->CandidateCells=PlaneCells*static_cast<uint64>(CellsPerAxis);
        OutStageTimings->ZBandRejectedCells=PlaneCells
            *static_cast<uint64>(CellsPerAxis-ActiveZ);
        if(bSkipProofBlockGrid)
        {
            for(int32 Z=GzLo;Z<=GzHi;++Z)
            for(int32 Y=-1;Y<=GridDim;++Y)
            for(int32 X=-1;X<=GridDim;++X)
            {
                const int32 Index=((Z+1)*MDim+(Y+1))*MDim+(X+1);
                if(RequiredGridSamples[Index]==0)
                    ++OutStageTimings->ProofSkippedDensitySamples;
            }
        }
    }
    const bool bCapturePreEditGrid = OutPreEditDensityGrid != nullptr
        && !bUseSharedDensityGrid
        && Step == 1 && CellsPerAxis == CHUNK_SIZE
        && CzLo == 0 && CzHi == CellsPerAxis - 1;
    if (bCapturePreEditGrid)
    {
        OutPreEditDensityGrid->OriginVoxels = OriginVoxels + FIntVector(-Step, -Step, -Step);
        OutPreEditDensityGrid->Step = Step;
        OutPreEditDensityGrid->Dim = MDim;
        OutPreEditDensityGrid->Samples.SetNumUninitialized(MDim * MDim * MDim);
    }
    // Bande de strate : seules les rangées Z réellement lues (cellules CzLo..CzHi + marges de
    // gradient) sont échantillonnées — le reste du tampon reste non initialisé et non lu.
    const uint64 StageDensityGridStartCycles=OutStageTimings?FPlatformTime::Cycles64():0;
    VoxelForgeCore::Linux77Profiler::FScopedStage l77GridFillStage(
        VoxelForgeCore::Linux77Profiler::EStage::GridFill);
    VoxelDensityProfile::FScopedTimer MesherDensityGridTimer(
        VoxelDensityProfile::EBucket::MesherDensityGrid);
#if WITH_EDITOR
    const bool bUseParameterizedBrickEvaluator =
        VoxelForgeRedefinedDensityField::IsParameterizedEnabled()
        && !bUseSharedDensityGrid && !bCapturePreEditGrid
        && OutCaptureGrid == nullptr;
#else
    constexpr bool bUseParameterizedBrickEvaluator = false;
#endif
    const bool bUseOperatorBlock = !bUseSharedDensityGrid
        && !bUseParameterizedBrickEvaluator
        && GVoxelForgeUseOperatorBlock != 0;
    bool bGridFillParallel=false;
    if (bUseParameterizedBrickEvaluator)
    {
#if WITH_EDITOR
        static std::atomic<bool> bLoggedParameterizedBrickRoute{false};
        if (!bLoggedParameterizedBrickRoute.exchange(true, std::memory_order_relaxed))
        {
            UE_LOG(LogTemp, Display,
                TEXT("[VoxelForgeResearchBrickEvaluator] mesher_grid=1 brick_samples=%d detail_scale=%.6f"),
                VoxelForgeRedefinedDensityField::ParameterizedBrickSamples,
                VoxelForgeRedefinedDensityField::GetDetailScale());
        }
        constexpr int32 ResearchBrickSamples =
            VoxelForgeRedefinedDensityField::ParameterizedBrickSamples;
        const float DetailScale = VoxelForgeRedefinedDensityField::GetDetailScale();
        TArray<float> DetailBrick;
        TArray<float> BrickValues;
        TArray<float> ColumnValues;
        TArray<uint8> ProofKinds;
        VoxelForgeRedefinedDensityField::FParameterizedBrickContext BrickContext;
        const int32 GridZSamples = GzHi - GzLo + 1;
        const int32 ProofZCount = FMath::DivideAndRoundUp(GridZSamples, ResearchBrickSamples);
        ProofKinds.SetNumUninitialized(ProofZCount);
        for (int32 BrickY = 0; BrickY < MDim; BrickY += ResearchBrickSamples)
        for (int32 BrickX = 0; BrickX < MDim; BrickX += ResearchBrickSamples)
        {
            if (ShouldAbortWork())
            {
                return FVoxelMeshData();
            }
            const int32 BrickSizeX = FMath::Min(ResearchBrickSamples, MDim - BrickX);
            const int32 BrickSizeY = FMath::Min(ResearchBrickSamples, MDim - BrickY);
            for (int32 ProofZ = 0; ProofZ < ProofZCount; ++ProofZ)
            {
                if (ShouldAbortWork())
                {
                    return FVoxelMeshData();
                }
                const int32 BrickZ = ProofZ * ResearchBrickSamples;
                const FIntVector ProofSize(
                    BrickSizeX, BrickSizeY,
                    FMath::Min(ResearchBrickSamples, GridZSamples - BrickZ));
                const FIntVector ProofOrigin = OriginVoxels + FIntVector(
                    (BrickX - 1) * Step,
                    (BrickY - 1) * Step,
                    (GzLo + BrickZ) * Step);

                bool bUniformAir = false;
                if (VoxelForgeRedefinedDensityField::TryClassifyParameterizedBrick(
                        *Generator, ProofOrigin, ProofSize, bUniformAir, DetailScale,
                        nullptr, nullptr, Step))
                {
                    ProofKinds[ProofZ] = bUniformAir ? 1 : 2;
                    const float UniformDensity = bUniformAir
                        ? VoxelForgeRedefinedDensityField::ParameterizedUniformDensity
                        : -VoxelForgeRedefinedDensityField::ParameterizedUniformDensity;
                    for (int32 Z = 0; Z < ProofSize.Z; ++Z)
                    for (int32 Y = 0; Y < ProofSize.Y; ++Y)
                    for (int32 X = 0; X < ProofSize.X; ++X)
                    {
                        const int32 GX = BrickX + X - 1;
                        const int32 GY = BrickY + Y - 1;
                        const int32 GZ = GzLo + BrickZ + Z;
                        const int32 LocalIndex =
                            ((GZ + 1) * MDim + (GY + 1)) * MDim + (GX + 1);
                        DensityGrid[LocalIndex] = UniformDensity;
                    }
                }
                else
                {
                    ProofKinds[ProofZ] = 0;
                }
            }

            for (int32 ProofZ = 0; ProofZ < ProofZCount;)
            {
                if (ProofKinds[ProofZ] != 0)
                {
                    ++ProofZ;
                    continue;
                }
                if (ShouldAbortWork())
                {
                    return FVoxelMeshData();
                }

                const int32 FirstMixedProofZ = ProofZ;
                const int32 BrickZ = FirstMixedProofZ * ResearchBrickSamples;
                const int32 FirstWorldZ = OriginVoxels.Z + (GzLo + BrickZ) * Step;
                const int32 FirstChunkZ = FMath::FloorToInt(
                    static_cast<float>(FirstWorldZ) / CHUNK_SIZE);
                int32 RunProofCount = 1;
                int32 RunSampleCount = FMath::Min(
                    ResearchBrickSamples, GridZSamples - BrickZ);
                while (ProofZ + RunProofCount < ProofZCount
                    && ProofKinds[ProofZ + RunProofCount] == 0)
                {
                    const int32 CandidateSampleCount = FMath::Min(
                        (RunProofCount + 1) * ResearchBrickSamples,
                        GridZSamples - BrickZ);
                    const int32 CandidateEndZ = FirstWorldZ
                        + (CandidateSampleCount - 1) * Step;
                    const int32 CandidateChunkZ = FMath::FloorToInt(
                        static_cast<float>(CandidateEndZ) / CHUNK_SIZE);
                    if (CandidateSampleCount > FMath::Max(1, CHUNK_SIZE / Step)
                        || CandidateChunkZ != FirstChunkZ)
                    {
                        break;
                    }
                    ++RunProofCount;
                    RunSampleCount = CandidateSampleCount;
                }

                const FIntVector BrickSize(BrickSizeX, BrickSizeY, RunSampleCount);
                const FIntVector BrickOrigin = OriginVoxels + FIntVector(
                    (BrickX - 1) * Step,
                    (BrickY - 1) * Step,
                    (GzLo + BrickZ) * Step);
                VoxelForgeRedefinedDensityField::EvaluateParameterizedDensityBrick(
                    *Generator, BrickOrigin, BrickSize, DetailScale, BrickContext,
                    DetailBrick, ColumnValues, BrickValues, Step);
                for (int32 Z = 0; Z < BrickSize.Z; ++Z)
                for (int32 Y = 0; Y < BrickSize.Y; ++Y)
                for (int32 X = 0; X < BrickSize.X; ++X)
                {
                    const int32 GX = BrickX + X - 1;
                    const int32 GY = BrickY + Y - 1;
                    const int32 GZ = GzLo + BrickZ + Z;
                    const int32 LocalIndex =
                        ((GZ + 1) * MDim + (GY + 1)) * MDim + (GX + 1);
                    const float Density=BrickValues[(Z * BrickSize.Y + Y) * BrickSize.X + X];
                    DensityGrid[LocalIndex]=Density;
                }
                if (OutDensitySampleCount != nullptr)
                {
                    *OutDensitySampleCount += static_cast<int64>(
                        BrickSize.X * BrickSize.Y * BrickSize.Z);
                }
                ProofZ += RunProofCount;
            }
        }
#endif
    }
    else
    {
        auto BuildAxisGroups = [&](int32 AxisOrigin, int32 RegionSize)
        {
            std::vector<std::vector<int32>> Groups;
            int32 PreviousRegion = INT32_MIN;
            for (int32 Index = -1; Index <= GridDim; ++Index)
            {
                const int64 WorldCoordinate = static_cast<int64>(AxisOrigin)
                    + static_cast<int64>(Index) * Step;
                const int32 Region = FMath::FloorToInt(
                    static_cast<float>(WorldCoordinate) / static_cast<float>(RegionSize));
                if (Groups.empty() || Region != PreviousRegion)
                {
                    Groups.emplace_back();
                    PreviousRegion = Region;
                }
                Groups.back().push_back(Index);
            }
            return Groups;
        };
        const int32 RoomCacheRegionSize = CHUNK_SIZE
            * (VoxelGenLOD::IsUseFusedEvaluatorEnabled() ? 4 : 1);
        const auto XGroups = BuildAxisGroups(OriginVoxels.X,RoomCacheRegionSize);
        const auto YGroups = BuildAxisGroups(OriginVoxels.Y,RoomCacheRegionSize);
        auto SampleGridPoint = [&](int32 gx, int32 gy, int32 gz,
                                  FVoxelMarchingCubesStageTimings* StageTimings,
                                  int64* DensitySampleCount,float* GridData)
        {
            const int32 LocalIndex = ((gz + 1) * MDim + (gy + 1)) * MDim + (gx + 1);
            if(bSkipProofBlockGrid&&RequiredGridSamples[LocalIndex]==0)return;
            if (bUseSharedDensityGrid)
            {
                const int32 SharedX = SharedDeltaX + gx;
                const int32 SharedY = SharedDeltaY + gy;
                const int32 SharedZ = SharedDeltaZ + gz;
                GridData[LocalIndex] = ReuseGrid->Samples[
                    (SharedZ * ReuseGrid->Dim + SharedY) * ReuseGrid->Dim + SharedX];
            }
            else
            {
                const int32 ReachBlockX = FMath::Clamp(gx, 0, CellsPerAxis - 1)
                    / VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockY = FMath::Clamp(gy, 0, CellsPerAxis - 1)
                    / VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockZ = FMath::Clamp(gz, 0, CellsPerAxis - 1)
                    / VoxelGenLOD::TileReachBlockCells;
                VoxelGenLOD::TileReachBlockIndex =
                    ((ReachBlockZ * ReachBlocksPerAxis) + ReachBlockY)
                        * ReachBlocksPerAxis + ReachBlockX;
                float PreEditDensity = 0.0f;
                const int32 SampleX=OriginVoxels.X+gx*Step;
                const int32 SampleY=OriginVoxels.Y+gy*Step;
                const int32 SampleZ=OriginVoxels.Z+gz*Step;
                const bool bBoundarySample=gx<0||gx>GridDim||gy<0||gy>GridDim
                    ||gz<0||gz>GridDim;
                VoxelForgeCore::Linux77Profiler::FScopedHalo l77Halo(bBoundarySample);
                VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller islandCaller(
                    bBoundarySample
                        ?VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Boundary
                        :VoxelForgeCore::floating_island_detail::EIslandFieldCaller::GridFill);
                GridData[LocalIndex] = ShapeWorldContext
                    ?ShapeWorldContext->DensityAtLattice(SampleX,SampleY,SampleZ)
                    :Generator->GetDensityAtLattice(SampleX,SampleY,SampleZ,
                        bCapturePreEditGrid ? &PreEditDensity : nullptr);
                if(StageTimings)
                {
                    if(bBoundarySample)++StageTimings->DensityGridHaloSamples;
                    else ++StageTimings->DensityGridInteriorSamples;
                }
                if (bCapturePreEditGrid)
                {
                    OutPreEditDensityGrid->Samples[LocalIndex] = PreEditDensity;
                }
                if (DensitySampleCount)
                {
                    ++(*DensitySampleCount);
                }
            }
        };
        // Floating-island LOD0 points reuse the same island/XY column keys across vertical
        // samples; large coarse tiles also revisit the same XY window across tall Z runs.
        // SampleGridPoint writes by lattice index, so output bits do not depend on traversal order.
        const int64 TileExtentVoxels = static_cast<int64>(Step) * CellsPerAxis;
        const float TileMiddleZ=static_cast<float>(OriginVoxels.Z)
            +static_cast<float>(TileExtentVoxels)*0.5f;
        const bool bFloatingIslandLod0=Step==1&&ShapeWorldContext
            &&ShapeWorldContext->IsFloatingIslandFieldAtZ(TileMiddleZ);
        const bool bUseVerticalCacheTraversal=bFloatingIslandLod0
            ||TileExtentVoxels>CHUNK_SIZE*16;
        const FStandaloneShapeWorld* ActiveShapeWorld=ShapeWorldContext;
        const bool bUseShapeWorldZBatch=ActiveShapeWorld!=nullptr
            &&ActiveShapeWorld->IsAnalyticOnly()
            &&!bUseSharedDensityGrid&&!bCapturePreEditGrid&&!bSkipProofBlockGrid;
        const bool bUseShapeWorldXBatch=ActiveShapeWorld!=nullptr
            &&!bUseSharedDensityGrid&&!bCapturePreEditGrid
            &&(!bSkipProofBlockGrid||bStandaloneEnableActiveLaneBatching);
        const auto SampleShapeWorldXBatch=[&](int32 GX0,int32 BatchCount,int32 gy,int32 gz,
            FVoxelMarchingCubesStageTimings* StageTimings,int64* DensitySampleCount,
            float* GridData)
        {
            float Densities[8]{};
            const bool bBoundarySample=GX0<0||GX0+BatchCount-1>GridDim
                ||gy<0||gy>GridDim||gz<0||gz>GridDim;
            VoxelForgeCore::Linux77Profiler::FScopedHalo l77Halo(bBoundarySample);
            VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller islandCaller(
                bBoundarySample
                    ?VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Boundary
                    :VoxelForgeCore::floating_island_detail::EIslandFieldCaller::GridFill);
            ActiveShapeWorld->DensityAtLatticeX8(OriginVoxels.X+GX0*Step,Step,
                OriginVoxels.Y+gy*Step,OriginVoxels.Z+gz*Step,BatchCount,Densities);
            if(StageTimings)
            {
                ++StageTimings->X8DensityCalls;
                StageTimings->X8ActiveLanes+=static_cast<uint64>(BatchCount);
                for(int32 Lane=0;Lane<BatchCount;++Lane)
                {
                    const int32 SampleGX=GX0+Lane;
                    if(SampleGX<0||SampleGX>GridDim||gy<0||gy>GridDim
                        ||gz<0||gz>GridDim)
                        ++StageTimings->DensityGridHaloSamples;
                    else ++StageTimings->DensityGridInteriorSamples;
                }
            }
            for(int32 Lane=0;Lane<BatchCount;++Lane)
            {
                const int32 GX=GX0+Lane;
                const int32 LocalIndex=((gz+1)*MDim+(gy+1))*MDim+(GX+1);
                GridData[LocalIndex]=Densities[Lane];
            }
            const int32 LastGX=GX0+BatchCount-1;
            const int32 ReachBlockX=FMath::Clamp(LastGX,0,CellsPerAxis-1)
                /VoxelGenLOD::TileReachBlockCells;
            const int32 ReachBlockY=FMath::Clamp(gy,0,CellsPerAxis-1)
                /VoxelGenLOD::TileReachBlockCells;
            const int32 ReachBlockZ=FMath::Clamp(gz,0,CellsPerAxis-1)
                /VoxelGenLOD::TileReachBlockCells;
            VoxelGenLOD::TileReachBlockIndex=
                ((ReachBlockZ*ReachBlocksPerAxis)+ReachBlockY)*ReachBlocksPerAxis+ReachBlockX;
            if(DensitySampleCount)*DensitySampleCount+=BatchCount;
        };
        const int32 GridZSampleCount=GzHi-GzLo+1;
        const int32 GridYSampleCount=GridDim+2;
        const bool bParallelizeGridRows=bUseShapeWorldZBatch
            ||bUseVerticalCacheTraversal||!bUseShapeWorldXBatch;
        const std::size_t GridFillLaneCount=std::min<std::size_t>({
            GridFillExecution.WorkerLimit,16,
            static_cast<std::size_t>(FMath::Max(0,
                bParallelizeGridRows?GridYSampleCount:GridZSampleCount))});
        const bool bUseParallelShapeWorldGrid=GridFillExecution.ParallelFor!=nullptr
            &&GridFillLaneCount>1&&GridFillExecution.BeginWorkerContext!=nullptr
            &&GridFillExecution.EndWorkerContext!=nullptr
            &&!bUseSharedDensityGrid&&!bCapturePreEditGrid&&!bSkipProofBlockGrid
            &&OutCaptureGrid==nullptr;
        if(bUseOperatorBlock&&!bUseParallelShapeWorldGrid)
        {
            // The grid includes one sample of halo on every side. The prepared operator block
            // is thread-local, so a parallel fill below opens one block per worker slab.
            Generator->BeginDensityBlock(
                OriginVoxels+FIntVector(-Step,-Step,GzLo*Step),
                Step,MDim,MDim,GzHi-GzLo+1);
        }
        if(bUseParallelShapeWorldGrid)
        {
            bGridFillParallel=true;
            const FVoxelMesherThreadState TileThreadState=FVoxelMesherThreadState::Capture();
            std::array<int64,16> LaneSampleCounts{};
            std::array<FVoxelMarchingCubesStageTimings,16> LaneTimings{};
            std::atomic<bool> bAborted{false};
            auto RunLane=[&](std::size_t Lane,std::size_t Worker)
            {
                void* WorkerState=GridFillExecution.BeginWorkerContext(
                    GridFillExecution.WorkerContext,Worker);
                int64& LaneSamples=LaneSampleCounts[Lane];
                FVoxelMarchingCubesStageTimings& StageTimings=LaneTimings[Lane];
                {
                    FVoxelMesherThreadStateScope ThreadStateScope(TileThreadState);
                    FStandaloneShapeWorld::FScopedPassageDensityQueryCaller queryCaller(
                        EPassageDensityQueryCaller::DensityGrid);
                    if(bUseOperatorBlock)
                        Generator->BeginDensityBlock(
                            OriginVoxels+FIntVector(-Step,-Step,GzLo*Step),Step,
                            MDim,MDim,GridZSampleCount);
                    if(bUseShapeWorldZBatch)
                    {
                        const int32 FirstGy=-1+static_cast<int32>(
                            (static_cast<std::size_t>(GridYSampleCount)*Lane)/GridFillLaneCount);
                        const int32 EndGy=-1+static_cast<int32>(
                            (static_cast<std::size_t>(GridYSampleCount)*(Lane+1))/GridFillLaneCount);
                        for(int32 gy=FirstGy;gy<EndGy;++gy)
                        {
                            if(ShouldAbortWork())
                            {bAborted.store(true,std::memory_order_relaxed);break;}
                            for(const std::vector<int32>& XGroup:XGroups)
                            for(const int32 gx:XGroup)
                            {
                                if(ShouldAbortWork())
                                {
                                    bAborted.store(true,std::memory_order_relaxed);
                                    break;
                                }
                                float Densities[CHUNK_SIZE+3]{};
                                const bool bBoundarySample=gx<0||gx>GridDim
                                    ||gy<0||gy>GridDim||GzLo<0||GzHi>GridDim;
                                VoxelForgeCore::Linux77Profiler::FScopedHalo l77Halo(
                                    bBoundarySample);
                                VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller
                                    islandCaller(bBoundarySample
                                        ?VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Boundary
                                        :VoxelForgeCore::floating_island_detail::EIslandFieldCaller::GridFill);
                                ActiveShapeWorld->DensityAtLatticeZRun(
                                    OriginVoxels.X+gx*Step,OriginVoxels.Y+gy*Step,
                                    OriginVoxels.Z+GzLo*Step,Step,GridZSampleCount,Densities);
                                for(int32 ZLane=0;ZLane<GridZSampleCount;++ZLane)
                                {
                                    const int32 gz=GzLo+ZLane;
                                    const int32 LocalIndex=
                                        ((gz+1)*MDim+(gy+1))*MDim+(gx+1);
                                    DensityGridData[LocalIndex]=Densities[ZLane];
                                    if(gx<0||gx>GridDim||gy<0||gy>GridDim
                                        ||gz<0||gz>GridDim)
                                        ++StageTimings.DensityGridHaloSamples;
                                    else ++StageTimings.DensityGridInteriorSamples;
                                }
                                LaneSamples+=GridZSampleCount;
                                const int32 ReachBlockX=FMath::Clamp(gx,0,CellsPerAxis-1)
                                    /VoxelGenLOD::TileReachBlockCells;
                                const int32 ReachBlockY=FMath::Clamp(gy,0,CellsPerAxis-1)
                                    /VoxelGenLOD::TileReachBlockCells;
                                const int32 ReachBlockZ=FMath::Clamp(GzHi,0,CellsPerAxis-1)
                                    /VoxelGenLOD::TileReachBlockCells;
                                VoxelGenLOD::TileReachBlockIndex=
                                    ((ReachBlockZ*ReachBlocksPerAxis)+ReachBlockY)
                                        *ReachBlocksPerAxis+ReachBlockX;
                            }
                            if(bAborted.load(std::memory_order_relaxed))break;
                        }
                    }
                    else if(bUseVerticalCacheTraversal||!bUseShapeWorldXBatch)
                    {
                        const int32 FirstGy=-1+static_cast<int32>(
                            (static_cast<std::size_t>(GridYSampleCount)*Lane)/GridFillLaneCount);
                        const int32 EndGy=-1+static_cast<int32>(
                            (static_cast<std::size_t>(GridYSampleCount)*(Lane+1))/GridFillLaneCount);
                        for(int32 gy=FirstGy;gy<EndGy;++gy)
                        {
                            if(ShouldAbortWork())
                            {bAborted.store(true,std::memory_order_relaxed);break;}
                            for(const std::vector<int32>& XGroup:XGroups)
                            for(const int32 gx:XGroup)
                            {
                                for(int32 gz=GzLo;gz<=GzHi;++gz)
                                {
                                    if(ShouldAbortWork())
                                    {bAborted.store(true,std::memory_order_relaxed);break;}
                                    SampleGridPoint(gx,gy,gz,&StageTimings,&LaneSamples,
                                        DensityGridData);
                                }
                                if(bAborted.load(std::memory_order_relaxed))break;
                            }
                            if(bAborted.load(std::memory_order_relaxed))break;
                        }
                    }
                    else
                    {
                        const int32 FirstGz=GzLo+static_cast<int32>(
                            (static_cast<std::size_t>(GridZSampleCount)*Lane)/GridFillLaneCount);
                        const int32 EndGz=GzLo+static_cast<int32>(
                            (static_cast<std::size_t>(GridZSampleCount)*(Lane+1))/GridFillLaneCount);
                        for(int32 gz=FirstGz;gz<EndGz;++gz)
                        {
                            if(ShouldAbortWork())
                            {bAborted.store(true,std::memory_order_relaxed);break;}
                            for(const std::vector<int32>& YGroup:YGroups)
                            for(const std::vector<int32>& XGroup:XGroups)
                            for(const int32 gy:YGroup)
                            {
                                if(ShouldAbortWork())
                                {bAborted.store(true,std::memory_order_relaxed);break;}
                                for(std::size_t Offset=0;Offset<XGroup.size();)
                                {
                                    const int32 BatchCount=static_cast<int32>(std::min<std::size_t>(
                                        8,XGroup.size()-Offset));
                                    const int32 GX0=XGroup[Offset];
                                    bool bContiguous=true;
                                    for(int32 BatchLane=1;BatchLane<BatchCount;++BatchLane)
                                        bContiguous&=XGroup[
                                            Offset+static_cast<std::size_t>(BatchLane)]==GX0+BatchLane;
                                    if(!bContiguous)
                                    {
                                        for(int32 BatchLane=0;BatchLane<BatchCount;++BatchLane)
                                            SampleGridPoint(XGroup[
                                                Offset+static_cast<std::size_t>(BatchLane)],gy,gz,
                                                &StageTimings,&LaneSamples,DensityGridData);
                                    }
                                    else
                                        SampleShapeWorldXBatch(GX0,BatchCount,gy,gz,
                                            &StageTimings,&LaneSamples,DensityGridData);
                                    Offset+=static_cast<std::size_t>(BatchCount);
                                }
                            }
                            if(bAborted.load(std::memory_order_relaxed))break;
                        }
                    }
                    if(bUseOperatorBlock)Generator->EndDensityBlock();
                }
                GridFillExecution.EndWorkerContext(
                    GridFillExecution.WorkerContext,WorkerState,Worker);
            };
            struct FDispatchContext{decltype(RunLane)* Run=nullptr;};
            FDispatchContext DispatchContext{&RunLane};
            const auto DispatchLane=[](void* RawContext,std::size_t Lane,std::size_t Worker)
            {
                auto& Context=*static_cast<FDispatchContext*>(RawContext);
                (*Context.Run)(Lane,Worker);
            };
            GridFillExecution.ParallelFor(GridFillExecution.SchedulerContext,
                GridFillLaneCount,GridFillLaneCount,&DispatchContext,DispatchLane);
            if(OutDensitySampleCount)
                for(std::size_t Lane=0;Lane<GridFillLaneCount;++Lane)
                    *OutDensitySampleCount+=LaneSampleCounts[Lane];
            if(OutStageTimings)
            {
                OutStageTimings->GridFillWorkerCount=static_cast<uint32>(GridFillLaneCount);
                for(std::size_t Lane=0;Lane<GridFillLaneCount;++Lane)
                {
                    OutStageTimings->DensityGridInteriorSamples+=
                        LaneTimings[Lane].DensityGridInteriorSamples;
                    OutStageTimings->DensityGridHaloSamples+=
                        LaneTimings[Lane].DensityGridHaloSamples;
                    OutStageTimings->X8DensityCalls+=LaneTimings[Lane].X8DensityCalls;
                    OutStageTimings->X8ActiveLanes+=LaneTimings[Lane].X8ActiveLanes;
                }
            }
            if(bAborted.load(std::memory_order_relaxed))return FVoxelMeshData();
        }
        else if(bUseShapeWorldZBatch&&!bUseVerticalCacheTraversal)
        {
            for(const std::vector<int32>& YGroup:YGroups)
            for(const std::vector<int32>& XGroup:XGroups)
            for(const int32 gy:YGroup)
            for(const int32 gx:XGroup)
            {
                if(ShouldAbortWork())
                {
                    if(bUseOperatorBlock){Generator->EndDensityBlock();}
                    return FVoxelMeshData();
                }
                const int32 BatchCount=GzHi-GzLo+1;
                const int32 SampleX=OriginVoxels.X+gx*Step;
                const int32 SampleY=OriginVoxels.Y+gy*Step;
                float Densities[CHUNK_SIZE+3]{};
                const bool bBoundarySample=gx<0||gx>GridDim||gy<0||gy>GridDim
                    ||GzLo<0||GzHi>GridDim;
                VoxelForgeCore::Linux77Profiler::FScopedHalo l77Halo(bBoundarySample);
                VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller islandCaller(
                    bBoundarySample
                        ?VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Boundary
                        :VoxelForgeCore::floating_island_detail::EIslandFieldCaller::GridFill);
                ActiveShapeWorld->DensityAtLatticeZRun(SampleX,SampleY,
                    OriginVoxels.Z+GzLo*Step,Step,BatchCount,Densities);
                if(OutStageTimings)
                {
                    for(int32 Lane=0;Lane<BatchCount;++Lane)
                    {
                        const int32 LocalGz=GzLo+Lane;
                        if(gx<0||gx>GridDim||gy<0||gy>GridDim
                            ||LocalGz<0||LocalGz>GridDim)
                            ++OutStageTimings->DensityGridHaloSamples;
                        else ++OutStageTimings->DensityGridInteriorSamples;
                    }
                }
                for(int32 Lane=0;Lane<BatchCount;++Lane)
                {
                    const int32 LocalGz=GzLo+Lane;
                    const int32 LocalIndex=
                        ((LocalGz+1)*MDim+(gy+1))*MDim+(gx+1);
                    DensityGrid[LocalIndex]=Densities[Lane];
                }
                const int32 ReachBlockX=FMath::Clamp(gx,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockY=FMath::Clamp(gy,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockZ=FMath::Clamp(GzHi,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                VoxelGenLOD::TileReachBlockIndex=
                    ((ReachBlockZ*ReachBlocksPerAxis)+ReachBlockY)
                        *ReachBlocksPerAxis+ReachBlockX;
                if(OutDensitySampleCount)*OutDensitySampleCount+=BatchCount;
            }
        }
        else if (bUseVerticalCacheTraversal)
        {
            for (const std::vector<int32>& YGroup : YGroups)
            for (const std::vector<int32>& XGroup : XGroups)
            for (const int32 gy : YGroup)
            for (const int32 gx : XGroup)
            {
                for (int32 gz = GzLo; gz <= GzHi; ++gz)
                {
                    if (ShouldAbortWork())
                    {
                        if (bUseOperatorBlock) { Generator->EndDensityBlock(); }
                        return FVoxelMeshData();
                    }
                    SampleGridPoint(gx,gy,gz,OutStageTimings,OutDensitySampleCount,
                        DensityGridData);
                }
            }
        }
        else
        {
            for (int32 gz = GzLo; gz <= GzHi; ++gz)
            {
                if (ShouldAbortWork())
                {
                    if (bUseOperatorBlock) { Generator->EndDensityBlock(); }
                    return FVoxelMeshData();
                }
                for (const std::vector<int32>& YGroup : YGroups)
                for (const std::vector<int32>& XGroup : XGroups)
                for (const int32 gy : YGroup)
                {
                    if (ShouldAbortWork())
                    {
                        if (bUseOperatorBlock) { Generator->EndDensityBlock(); }
                        return FVoxelMeshData();
                    }
                    if(bUseShapeWorldXBatch)
                    {
                        for(size_t Offset=0;Offset<XGroup.size();)
                        {
                            if(bSkipProofBlockGrid)
                            {
                                const int32 GX0=XGroup[Offset];
                                const int32 LocalIndex=((gz+1)*MDim+(gy+1))*MDim+(GX0+1);
                                if(RequiredGridSamples[LocalIndex]==0)
                                {++Offset;continue;}
                                int32 BatchCount=1;
                                while(BatchCount<8&&Offset+static_cast<size_t>(BatchCount)<XGroup.size()
                                    &&XGroup[Offset+static_cast<size_t>(BatchCount)]==GX0+BatchCount)
                                {
                                    const int32 NextGX=XGroup[Offset+static_cast<size_t>(BatchCount)];
                                    const int32 NextIndex=((gz+1)*MDim+(gy+1))*MDim+(NextGX+1);
                                    if(RequiredGridSamples[NextIndex]==0)break;
                                    ++BatchCount;
                                }
                                if(BatchCount<2)SampleGridPoint(
                                    GX0,gy,gz,OutStageTimings,OutDensitySampleCount,
                                    DensityGridData);
                                else SampleShapeWorldXBatch(GX0,BatchCount,gy,gz,
                                    OutStageTimings,OutDensitySampleCount,DensityGridData);
                                Offset+=static_cast<size_t>(BatchCount);
                                continue;
                            }
                            const int32 BatchCount=static_cast<int32>(std::min<size_t>(
                                8,XGroup.size()-Offset));
                            const int32 GX0=XGroup[Offset];
                            bool bContiguous=true;
                            for(int32 Lane=1;Lane<BatchCount;++Lane)
                                bContiguous&=XGroup[Offset+static_cast<size_t>(Lane)]==GX0+Lane;
                            if(!bContiguous)
                            {
                                for(int32 Lane=0;Lane<BatchCount;++Lane)
                                    SampleGridPoint(XGroup[Offset+static_cast<size_t>(Lane)],gy,gz,
                                        OutStageTimings,OutDensitySampleCount,DensityGridData);
                                Offset+=static_cast<size_t>(BatchCount);
                                continue;
                            }

                            SampleShapeWorldXBatch(GX0,BatchCount,gy,gz,
                                OutStageTimings,OutDensitySampleCount,DensityGridData);
                            Offset+=static_cast<size_t>(BatchCount);
                        }
                    }
                    else
                    {
                        for(const int32 gx:XGroup)SampleGridPoint(
                            gx,gy,gz,OutStageTimings,OutDensitySampleCount,DensityGridData);
                    }
                }
            }
        }
    }
    if (VoxelDensityProfile::AreCountersEnabled())
    {
        const uint64 GridSampleCount = static_cast<uint64>(GzHi - GzLo + 1)
            * static_cast<uint64>(MDim)
            * static_cast<uint64>(MDim);
        VoxelDensityProfile::AddCounter(
            VoxelDensityProfile::ECounter::MesherDensityGridBytes,
            GridSampleCount * sizeof(float));
        if (bUseSharedDensityGrid)
        {
            VoxelDensityProfile::AddCounter(
                VoxelDensityProfile::ECounter::MesherSharedGridReadBytes,
                GridSampleCount * sizeof(float));
        }
    }
    if (bUseOperatorBlock&&!bGridFillParallel)
    {
        Generator->EndDensityBlock();
    }
    if (OutDebugDensitySamples)
    {
        // Hash/debug runs need the complete reference density lattice. Let MC exercise the
        // sparse block path, then evaluate omitted samples before exporting that lattice.
        if(bSkipProofBlockGrid)
        {
            for(int32 gz=GzLo;gz<=GzHi;++gz)
            for(int32 gy=-1;gy<=GridDim;++gy)
            for(int32 gx=-1;gx<=GridDim;++gx)
            {
                const int32 LocalIndex=((gz+1)*MDim+(gy+1))*MDim+(gx+1);
                if(RequiredGridSamples[LocalIndex]!=0)continue;
                RequiredGridSamples[LocalIndex]=1;
                const int32 ReachBlockX=FMath::Clamp(gx,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockY=FMath::Clamp(gy,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                const int32 ReachBlockZ=FMath::Clamp(gz,0,CellsPerAxis-1)
                    /VoxelGenLOD::TileReachBlockCells;
                VoxelGenLOD::TileReachBlockIndex=
                    ((ReachBlockZ*ReachBlocksPerAxis)+ReachBlockY)*ReachBlocksPerAxis+ReachBlockX;
                const int32 SampleX=OriginVoxels.X+gx*Step;
                const int32 SampleY=OriginVoxels.Y+gy*Step;
                const int32 SampleZ=OriginVoxels.Z+gz*Step;
                const bool bBoundarySample=gx<0||gx>GridDim||gy<0||gy>GridDim;
                VoxelForgeCore::Linux77Profiler::FScopedHalo l77Halo(bBoundarySample);
                VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller islandCaller(
                    bBoundarySample
                        ?VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Boundary
                        :VoxelForgeCore::floating_island_detail::EIslandFieldCaller::GridFill);
                DensityGrid[LocalIndex]=ShapeWorldContext
                    ?ShapeWorldContext->DensityAtLattice(SampleX,SampleY,SampleZ)
                    :Generator->GetDensityAtLattice(SampleX,SampleY,SampleZ);
                if(OutDensitySampleCount)++(*OutDensitySampleCount);
            }
        }
        const size_t SampleCount = static_cast<size_t>(GzHi - GzLo + 1)
            * static_cast<size_t>(MDim) * static_cast<size_t>(MDim);
        OutDebugDensitySamples->reserve(SampleCount);
        for (int32 gz = GzLo; gz <= GzHi; ++gz)
        for (int32 gy = -1; gy <= GridDim; ++gy)
        for (int32 gx = -1; gx <= GridDim; ++gx)
        {
            const int32 LocalIndex = ((gz + 1) * MDim + (gy + 1)) * MDim + (gx + 1);
            OutDebugDensitySamples->push_back(DensityGrid[LocalIndex]);
        }
    }
    MesherDensityGridTimer.End();
    l77GridFillStage.Restore();
    VoxelForgeCore::Linux77Profiler::FScopedStage l77MarchingCubesStage(
        VoxelForgeCore::Linux77Profiler::EStage::MarchingCubes);
    if(OutStageTimings)
        OutStageTimings->DensityGridMilliseconds=static_cast<double>(
            FPlatformTime::Cycles64()-StageDensityGridStartCycles)
            *FPlatformTime::GetSecondsPerCycle64()*1000.0;

    if (bUseProofBlocks&&(!bSkipProofBlockGrid||OutDebugDensitySamples!=nullptr))
    {
        // Finalize every block candidate against the exact density samples that the marching-
        // cubes loop will consume.  This is intentionally redundant with ClassifyTile: it is a
        // cheap soundness fence around the proof/cache hand-off.  Unknown, zero, or opposite-sign
        // vertices clear the candidate, so an unproven block pays the ordinary cell path.
        const int32 BlockCells = ProofBlockCells;
        for (int32 BlockZ = 0; BlockZ < BlocksPerAxis; ++BlockZ)
        {
            const int32 BlockMinZ = BlockZ * BlockCells;
            const int32 BlockMaxZ = FMath::Min(CellsPerAxis - 1,
                                               BlockMinZ + BlockCells - 1);
            if (BlockMaxZ < CzLo || BlockMinZ > CzHi)
            {
                continue;
            }
            for (int32 BlockY = 0; BlockY < BlocksPerAxis; ++BlockY)
            for (int32 BlockX = 0; BlockX < BlocksPerAxis; ++BlockX)
            {
                const int32 BlockIndex =
                    ((BlockZ * BlocksPerAxis) + BlockY) * BlocksPerAxis + BlockX;
                const uint8 Candidate = BlockSkip[BlockIndex];
                if (Candidate == 0)
                {
                    continue;
                }

                const int32 BlockMaxX = FMath::Min(CellsPerAxis - 1,
                                                   BlockX * BlockCells + BlockCells - 1);
                const int32 BlockMaxY = FMath::Min(CellsPerAxis - 1,
                                                   BlockY * BlockCells + BlockCells - 1);
                bool bUniform = true;
                for (int32 Z = BlockMinZ; Z <= BlockMaxZ + 1 && bUniform; ++Z)
                for (int32 Y = BlockY * BlockCells; Y <= BlockMaxY + 1 && bUniform; ++Y)
                for (int32 X = BlockX * BlockCells; X <= BlockMaxX + 1; ++X)
                {
                    const float Density = DensityGrid[
                        ((Z + 1) * MDim + (Y + 1)) * MDim + (X + 1)];
                    if (!VoxelMath::IsFinite(Density)
                        || Density == 0.0f
                        || (Candidate == 1 ? Density >= 0.0f : Density <= 0.0f))
                    {
                        bUniform = false;
                        break;
                    }
                }
                if (!bUniform)
                {
                    BlockSkip[BlockIndex] = 0;
                    continue;
                }

                if (VoxelDensityProfile::AreCountersEnabled())
                {
                    VoxelDensityProfile::AddCounter(
                        Candidate == 1
                            ? VoxelDensityProfile::ECounter::MesherBlockAllSolid
                            : VoxelDensityProfile::ECounter::MesherBlockAllAir);
                }
            }
        }
    }
    if(OutStageTimings&&bUseProofBlocks)
    {
        for(int32 BlockZ=0;BlockZ<BlocksPerAxis;++BlockZ)
        {
            const int32 CellMinZ=BlockZ*ProofBlockCells;
            const int32 CellMaxZ=FMath::Min(CellsPerAxis-1,CellMinZ+ProofBlockCells-1);
            const int32 ZCount=FMath::Max(0,FMath::Min(CzHi,CellMaxZ)
                -FMath::Max(CzLo,CellMinZ)+1);
            if(ZCount==0)continue;
            for(int32 BlockY=0;BlockY<BlocksPerAxis;++BlockY)
            for(int32 BlockX=0;BlockX<BlocksPerAxis;++BlockX)
            {
                const int32 Index=(BlockZ*BlocksPerAxis+BlockY)*BlocksPerAxis+BlockX;
                const uint8 Reason=BlockSkip[Index];
                if(Reason==0)continue;
                const int32 XCount=FMath::Min(CellsPerAxis-1,
                    (BlockX+1)*ProofBlockCells-1)-BlockX*ProofBlockCells+1;
                const int32 YCount=FMath::Min(CellsPerAxis-1,
                    (BlockY+1)*ProofBlockCells-1)-BlockY*ProofBlockCells+1;
                const uint64 Skipped=static_cast<uint64>(XCount)*YCount*ZCount;
                if(Reason==1)OutStageTimings->ProofSkippedSolidCells+=Skipped;
                else if(Reason==2)OutStageTimings->ProofSkippedAirCells+=Skipped;
            }
        }
    }

    // ── CAPTURE-DURING-MESHING ──
    // Si demandé et que la tuile est pleine résolution (CellsPerAxis==CHUNK_SIZE ⇒ Step==1<<Level,
    // donc chaque point de grille = exactement une cellule du clipmap de densité), on recopie les
    // CHUNK_SIZE³ points INTÉRIEURS (g=0..CHUNK_SIZE-1, on exclut le point frontière +1 — il
    // appartient à la tuile voisine — et l'anneau de marge ±1) dans OutCaptureGrid, quantifiés.
    // UVoxelDensityVolume réutilise ces octets au lieu de re-sampler GetDensityAt. Pure lecture de
    // DensityGrid : la forme de grille, la boucle deux passes, l'anneau de marge et la réutilisation
    // thread_local restent intacts (§8.10). Ordre X→Y→Z (x rapide) = layout attendu par l'ingest.
    if (OutCaptureGrid && CellsPerAxis == CHUNK_SIZE)
    {
        OutCaptureGrid->SetNumUninitialized(CHUNK_SIZE * CHUNK_SIZE * CHUNK_SIZE);
        uint8* Cap = OutCaptureGrid->GetData();
        int32 ci = 0;
        for (int32 gz = 0; gz < CHUNK_SIZE; ++gz)
        for (int32 gy = 0; gy < CHUNK_SIZE; ++gy)
        for (int32 gx = 0; gx < CHUNK_SIZE; ++gx)
        {
            Cap[ci++] = VF_QuantizeDensity(DensityGrid[((gz + 1) * MDim + (gy + 1)) * MDim + (gx + 1)]);
        }
    }

    if (OutPreEditDensityGrid && !bCapturePreEditGrid)
    {
        OutPreEditDensityGrid->Reset();
    }

    // Lecture grille (avec offset de marge) + gradient central depuis la grille.
    auto SampleG = [&](int32 gx, int32 gy, int32 gz) -> float
    {
        return DensityGrid[((gz + 1) * MDim + (gy + 1)) * MDim + (gx + 1)];
    };
    // Reuse exact grid-corner gradients between incident cells within this tile.
    // Reutiliser les gradients exacts des coins entre cellules voisines de cette tuile.
    static thread_local TArray<FVector> GradientMemo;
    static thread_local TArray<uint8> GradientMemoValid;
    const bool bUseGradientMemo=bStandaloneEnableGradientMemoization
        &&ShapeWorldContext!=nullptr&&Step==1;
    if(bUseGradientMemo)
    {
        GradientMemo.SetNumUninitialized(DensityGrid.Num());
        GradientMemoValid.SetNumUninitialized(DensityGrid.Num());
        // Clear every tile even when the reused TArray keeps the same length.
        // Effacer chaque tuile meme si le TArray reutilise conserve la meme taille.
        if(GradientMemoValid.Num()>0)
            FMemory::Memzero(GradientMemoValid.GetData(),
                static_cast<SIZE_T>(GradientMemoValid.Num())*sizeof(uint8));
    }
    auto SampleGradientG = [&](int32 gx, int32 gy, int32 gz) -> float
    {
        const float Cached = SampleG(gx, gy, gz);
        if (!ShapeWorldContext || Step != 1
            || !ShapeWorldContext->IsProxyBandFallback(Cached))
        {
            return Cached;
        }
        if(OutStageTimings)++OutStageTimings->GradientFallbackDensityCalls;
        FStandaloneShapeWorld::FScopedPassageDensityQueryCaller queryCaller(
            EPassageDensityQueryCaller::GradientNormalFallback);
        return ShapeWorldContext->DensityAtUnclamped(
            static_cast<float>(OriginVoxels.X + gx * Step),
            static_cast<float>(OriginVoxels.Y + gy * Step),
            static_cast<float>(OriginVoxels.Z + gz * Step));
    };
    auto GradAt = [&](int32 gx, int32 gy, int32 gz) -> FVector
    {
        if(OutStageTimings)++OutStageTimings->GradientCalls;
        const int32 CacheIndex=((gz+1)*MDim+(gy+1))*MDim+(gx+1);
        if(bUseGradientMemo&&GradientMemoValid[CacheIndex]!=0)
        {
            if(OutStageTimings)++OutStageTimings->GradientMemoHits;
            return GradientMemo[CacheIndex];
        }
        // Densité négative=solide, positive=air → le gradient pointe vers l'air (sortant).
        const FVector Gradient(
            SampleGradientG(gx + 1, gy, gz) - SampleGradientG(gx - 1, gy, gz),
            SampleGradientG(gx, gy + 1, gz) - SampleGradientG(gx, gy - 1, gz),
            SampleGradientG(gx, gy, gz + 1) - SampleGradientG(gx, gy, gz - 1));
        if(bUseGradientMemo)
        {
            GradientMemo[CacheIndex]=Gradient;
            GradientMemoValid[CacheIndex]=1;
        }
        return Gradient;
    };

    // DC is limited to tagged LOD0 patches. Coarse tiles and a connected patch that reaches a
    // tile edge use MC there, preserving the existing MC-to-MC clipmap transition and skirts.
    // Le DC reste limite aux surfaces LOD0 balisees; les patches connectes au bord restent MC.
    VoxelHybridSharpMesher::FResult HybridResult;
    const bool bHybridLod0Interior = Step == 1 && CellsPerAxis == CHUNK_SIZE
        && CzLo == 0 && CzHi == CellsPerAxis - 1 && !HybridSharpFeatures.IsEmpty();
    if (bUseHybridSharpMesher && bHybridLod0Interior)
    {
        HybridResult = VoxelHybridSharpMesher::Generate(
            CellsPerAxis, Step, IsoLevel, GeometryVoxelSize, OriginVoxels, ChunkWorldPos,
            HybridSharpFeatures,
            [&](int32 X, int32 Y, int32 Z) { return SampleG(X, Y, Z); },
            [&](int32 X, int32 Y, int32 Z) { return GradAt(X, Y, Z); });
    }
    const bool bHybridSharpSnap = bUseHybridSharpFeatureSnap && bHybridLod0Interior;
    const std::vector<VoxelHybridSharpMesher::FPreparedFeatureBox> PreparedSharpFeatures =
        bHybridSharpSnap
            ? VoxelHybridSharpMesher::PrepareFeatureBoxes(HybridSharpFeatures)
            : std::vector<VoxelHybridSharpMesher::FPreparedFeatureBox>{};

    //=========================================================================
    // ITÉRATION SUR LES CELLULES
    //=========================================================================
    // On itère sur les CELLULES (indices de grille), pas sur les voxels monde.
    // Coin i de la cellule = point de grille (cx+ox, cy+oy, cz+oz) ; coord voxel
    // monde = ce point × Step. LOD0 Step=1 → full res ; LOD1 Step=2 → ~4× moins.
    // F17 — les triangles sont émis dans DEUX seaux (sol / sky-cap, vote majoritaire des
    // classes de vertex) puis concaténés sol-puis-cap : RMC exige un run d'indices contigu
    // par polygroup. Géométrie inchangée au bit près — seul l'ORDRE des triangles bouge.
    static thread_local TArray<int32> GroundTris;
    static thread_local TArray<int32> CapTris;
    GroundTris.Reset();
    CapTris.Reset();
    const bool bProfileMesher = VoxelDensityProfile::IsCycleTimingEnabled();
    const bool bAuditMesherStages=OutStageTimings!=nullptr;
    uint64 CellClassificationCalls = 0;
    uint64 SurfaceCellCalls = 0;
    uint64 AuditCellClassificationCycles=0,AuditGradientCycles=0;
    uint64 AuditInterpolationCycles=0,AuditDedupStreamCycles=0;
    uint64 AuditGradientSamples=0,AuditInterpolationSamples=0,AuditDedupStreamSamples=0;
    struct FMidpointDensityKey
    {
        uint64 X=0,Y=0,Z=0;
        bool operator==(const FMidpointDensityKey&)const=default;
    };
    struct FMidpointDensityKeyHash
    {
        size_t operator()(const FMidpointDensityKey& Key)const noexcept
        {
            uint64 Hash=Key.X^std::rotl(Key.Y,21)^std::rotl(Key.Z,42);
            Hash^=Hash>>30;Hash*=0xbf58476d1ce4e5b9ull;
            Hash^=Hash>>27;Hash*=0x94d049bb133111ebull;
            Hash^=Hash>>31;
            return static_cast<size_t>(Hash);
        }
    };
    std::unordered_map<FMidpointDensityKey,float,FMidpointDensityKeyHash> MidpointDensityMemo;
    if(bStandaloneEnableExactMidpointDensityMemo&&ShapeWorldContext&&Step==1)
        MidpointDensityMemo.reserve(1024);
    uint8 CaseIndexRow[CHUNK_SIZE]{};
    uint8 SurfaceRow[CHUNK_SIZE]{};
    uint8 SampleClassificationRow[CHUNK_SIZE]{};
    for (int32 cz = CzLo; cz <= CzHi; cz++)            // bande de strate : cf. CzLo/CzHi plus haut
    {
        if (ShouldAbortWork()) return FVoxelMeshData();
        for (int32 cy = 0; cy < CellsPerAxis; cy++)
        {
            if (ShouldAbortWork()) return FVoxelMeshData();
            FMemory::Memzero(SampleClassificationRow, sizeof(SampleClassificationRow));
            const uint64 ClassificationRowStart = (bProfileMesher||bAuditMesherStages)
                ? FPlatformTime::Cycles64() : 0;
            int32 ClassificationRowCalls = 0;
            for (int32 GroupX = 0; GroupX < CellsPerAxis; GroupX += 8)
            {
                const int32 GroupCount = FMath::Min(8, CellsPerAxis - GroupX);
                bool bAnySkipped=false,bAllSkipped=bUseProofBlocks;
                if(bUseProofBlocks)
                {
                    for(int32 Lane=0;Lane<GroupCount;++Lane)
                    {
                        const int32 CellX=GroupX+Lane;
                        const int32 BlockIndex=((cz/ProofBlockCells)*BlocksPerAxis
                            +(cy/ProofBlockCells))*BlocksPerAxis+(CellX/ProofBlockCells);
                        const bool bSkipped=BlockSkip[BlockIndex]!=0;
                        bAnySkipped|=bSkipped;bAllSkipped&=bSkipped;
                    }
                }
                if(bAllSkipped)continue;

                if (GroupCount == 8&&!bAnySkipped)
                {
                    const float* Row00 = DensityGrid.GetData()
                        + ((cz + 1) * MDim + (cy + 1)) * MDim + GroupX + 1;
                    const float* Row01 = Row00 + MDim;
                    const float* Row10 = Row00 + MDim * MDim;
                    const float* Row11 = Row10 + MDim;
                    const __m256 Iso = _mm256_set1_ps(IsoLevel);
                    const int32 Masks[8] = {
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row00), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row00 + 1), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row01 + 1), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row01), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row10), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row10 + 1), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row11 + 1), Iso, _CMP_GE_OQ)),
                        _mm256_movemask_ps(_mm256_cmp_ps(_mm256_loadu_ps(Row11), Iso, _CMP_GE_OQ)),
                    };
                    for (int32 Lane = 0; Lane < 8; ++Lane)
                    {
                        const int32 CellX = GroupX + Lane;
                        int32 CaseIndex = 0;
                        for (int32 Corner = 0; Corner < 8; ++Corner)
                        {
                            CaseIndex |= ((Masks[Corner] >> Lane) & 1) << Corner;
                        }
                        CaseIndexRow[CellX] = static_cast<uint8>(CaseIndex);
                        SurfaceRow[CellX] = MCEdgeTable[CaseIndex] != 0;
                        SampleClassificationRow[CellX] =
                            (bProfileMesher&&VoxelDensityProfile::ShouldSample(
                                VoxelDensityProfile::EBucket::MesherCellClassification))
                            ||(bAuditMesherStages&&(CellClassificationCalls%8u)==0u);
                        ++ClassificationRowCalls;
                        ++CellClassificationCalls;
                    }
                }
                else
                {
                    for (int32 Lane = 0; Lane < GroupCount; ++Lane)
                    {
                        const int32 CellX = GroupX + Lane;
                        if(bUseProofBlocks
                            &&BlockSkip[((cz/ProofBlockCells)*BlocksPerAxis
                                +(cy/ProofBlockCells))*BlocksPerAxis
                                +(CellX/ProofBlockCells)]!=0)continue;
                        int32 CaseIndex = 0;
                        for (int32 Corner = 0; Corner < 8; ++Corner)
                        {
                            const float Density = SampleG(
                                CellX + CornerOffsets[Corner].X,
                                cy + CornerOffsets[Corner].Y,
                                cz + CornerOffsets[Corner].Z);
                            if (Density >= IsoLevel) CaseIndex |= 1 << Corner;
                        }
                        CaseIndexRow[CellX] = static_cast<uint8>(CaseIndex);
                        SurfaceRow[CellX] = MCEdgeTable[CaseIndex] != 0;
                        SampleClassificationRow[CellX] =
                            (bProfileMesher&&VoxelDensityProfile::ShouldSample(
                                VoxelDensityProfile::EBucket::MesherCellClassification))
                            ||(bAuditMesherStages&&(CellClassificationCalls%8u)==0u);
                        ++ClassificationRowCalls;
                        ++CellClassificationCalls;
                    }
                }
            }
            if ((bProfileMesher||bAuditMesherStages) && ClassificationRowCalls > 0)
            {
                const uint64 RowCycles = FPlatformTime::Cycles64() - ClassificationRowStart;
                if(bAuditMesherStages)AuditCellClassificationCycles+=RowCycles;
                if(bProfileMesher)
                {
                    const uint64 PerCellCycles = RowCycles
                        / static_cast<uint64>(ClassificationRowCalls);
                    for (int32 cx = 0; cx < CellsPerAxis; ++cx)
                    {
                        if (SampleClassificationRow[cx])
                        {
                            VoxelDensityProfile::AddSampledMeasurement(
                                VoxelDensityProfile::EBucket::MesherCellClassification,
                                PerCellCycles);
                        }
                    }
                }
            }
            for (int32 cx = 0; cx < CellsPerAxis; cx++)
            {
                if (bUseProofBlocks
                    && BlockSkip[((cz / ProofBlockCells) * BlocksPerAxis
                        + (cy / ProofBlockCells)) * BlocksPerAxis
                        + (cx / ProofBlockCells)] != 0)
                {
                    continue;
                }
                const bool bSampleCell = SampleClassificationRow[cx] != 0;
                float Densities[8];
                const int32 CaseIndex = CaseIndexRow[cx];
                const bool bHasSurface = SurfaceRow[cx] != 0;
                if (!bHasSurface) continue;  // Pas de surface ici → skip
                ++SurfaceCellCalls;
                // Density values are needed only after classification, for edge interpolation.
                for (int32 i = 0; i < 8; i++)
                {
                    Densities[i] = SampleG(cx + CornerOffsets[i].X,
                                           cy + CornerOffsets[i].Y,
                                           cz + CornerOffsets[i].Z);
                }

                // PASSE 2 : positions + gradients aux 8 coins (uniquement si surface présente).
                const uint64 GradientStart = bSampleCell
                    ? FPlatformTime::Cycles64() : 0;
                FVector Positions[8];
                FVector Gradients[8];
                for (int32 i = 0; i < 8; i++)
                {
                    const int32 GX = cx + CornerOffsets[i].X;
                    const int32 GY = cy + CornerOffsets[i].Y;
                    const int32 GZ = cz + CornerOffsets[i].Z;
                    Positions[i] = ChunkWorldPos
                        + FVector(GX * Step, GY * Step, GZ * Step) * GeometryVoxelSize;
                    Gradients[i] = GradAt(GX, GY, GZ);
                }
                if (bSampleCell)
                {
                    const uint64 StageCycles=FPlatformTime::Cycles64()-GradientStart;
                    if(bProfileMesher)VoxelDensityProfile::AddSampledMeasurement(
                        VoxelDensityProfile::EBucket::MesherGradientNormals,StageCycles);
                    if(bAuditMesherStages)
                    {AuditGradientCycles+=StageCycles;++AuditGradientSamples;}
                }

                // Interpolation des positions + normales sur les arêtes traversées. t = point de
                // traversée de l'iso entre les deux coins (clampé, milieu si densités quasi-égales) ;
                // la normale interpole les gradients de coin par le même t.
                const uint64 VertexInterpolationStart = bSampleCell
                    ? FPlatformTime::Cycles64() : 0;
                FVector EdgeVertices[12];
                FVector EdgeNormals[12];
                uint64 EdgeSplitKeys[12]{};
                for (int32 i = 0; i < 12; i++)
                {
                    if (MCEdgeTable[CaseIndex] & (1 << i))
                    {
                        const int32 A = EdgeCorners[i][0];
                        const int32 B = EdgeCorners[i][1];
                        const float D1 = Densities[A], D2 = Densities[B];
                        float T = (FMath::Abs(D2 - D1) < KINDA_SMALL_NUMBER)
                            ? 0.5f
                            : FMath::Clamp((IsoLevel - D1) / (D2 - D1), 0.0f, 1.0f);
                        if(ShapeWorldContext&&Step==1&&T>=0.375f&&T<=0.625f)
                        {
                            const FVector Midpoint=(Positions[A]+Positions[B])*0.5f;
                            if(OutStageTimings)++OutStageTimings->EdgeMidpointDensityCalls;
                            const double MidX=Midpoint.X/GeometryVoxelSize;
                            const double MidY=Midpoint.Y/GeometryVoxelSize;
                            const double MidZ=Midpoint.Z/GeometryVoxelSize;
                            const auto EvaluateExactMidpoint=[&]()
                            {
                                FStandaloneShapeWorld::FScopedPassageDensityQueryCaller
                                    queryCaller(EPassageDensityQueryCaller::MidpointCorrection);
                                return ShapeWorldContext->DensityAtUnclamped(
                                    MidX,MidY,MidZ);
                            };
                            float MidDensity=0.0f;
                            if(bStandaloneEnableExactMidpointDensityMemo)
                            {
                                const FMidpointDensityKey Key{std::bit_cast<uint64>(MidX),
                                    std::bit_cast<uint64>(MidY),std::bit_cast<uint64>(MidZ)};
                                const auto Found=MidpointDensityMemo.find(Key);
                                if(Found!=MidpointDensityMemo.end())
                                {
                                    MidDensity=Found->second;
                                    if(OutStageTimings)
                                        ++OutStageTimings->EdgeMidpointDensityMemoHits;
                                }
                                else
                                {
                                    MidDensity=EvaluateExactMidpoint();
                                    MidpointDensityMemo.emplace(Key,MidDensity);
                                }
                            }
                            else MidDensity=EvaluateExactMidpoint();
                            if(FMath::Abs(MidDensity-IsoLevel)<1.0e-6f)
                                T=0.5f;
                            else if((MidDensity<IsoLevel)==(D1<IsoLevel))
                            {
                                const float Denominator=D2-MidDensity;
                                if(FMath::Abs(Denominator)>KINDA_SMALL_NUMBER)
                                    T=0.5f+0.5f*FMath::Clamp(
                                        (IsoLevel-MidDensity)/Denominator,0.0f,1.0f);
                            }
                            else
                            {
                                const float Denominator=MidDensity-D1;
                                if(FMath::Abs(Denominator)>KINDA_SMALL_NUMBER)
                                    T=0.5f*FMath::Clamp(
                                        (IsoLevel-D1)/Denominator,0.0f,1.0f);
                            }
                        }
                        EdgeVertices[i] = Positions[A] + T * (Positions[B] - Positions[A]);
                        EdgeNormals[i]  = Gradients[A] + T * (Gradients[B] - Gradients[A]);
                        if (bHybridSharpSnap)
                        {
                            const int32 Axis = CornerOffsets[A].X != CornerOffsets[B].X ? 0
                                : (CornerOffsets[A].Y != CornerOffsets[B].Y ? 1 : 2);
                            const FIntVector LocalStart(cx + CornerOffsets[A].X,
                                cy + CornerOffsets[A].Y, cz + CornerOffsets[A].Z);
                            // All cells incident to this lattice edge must stay interior so
                            // the shared root cannot disagree across a tile-border cell.
                            const bool bCanSnapEdge = bHybridSharpSnap
                                && LocalStart.X >= 2 && LocalStart.Y >= 2 && LocalStart.Z >= 2
                                && LocalStart.X <= CellsPerAxis - 2
                                && LocalStart.Y <= CellsPerAxis - 2
                                && LocalStart.Z <= CellsPerAxis - 2;
                            FVector AnalyticNormal;
                            if (bCanSnapEdge && VoxelHybridSharpMesher::TrySnapFeatureEdge(
                                    PreparedSharpFeatures, OriginVoxels, LocalStart, Axis,
                                    ChunkWorldPos, GeometryVoxelSize, EdgeVertices[i],
                                    AnalyticNormal, EdgeSplitKeys[i]))
                            {
                                EdgeNormals[i] = AnalyticNormal;
                            }
                        }
                    }
                }
                if (bSampleCell)
                {
                    const uint64 StageCycles=FPlatformTime::Cycles64()-VertexInterpolationStart;
                    if(bProfileMesher)VoxelDensityProfile::AddSampledMeasurement(
                        VoxelDensityProfile::EBucket::MesherVertexInterpolation,StageCycles);
                    if(bAuditMesherStages)
                    {AuditInterpolationCycles+=StageCycles;++AuditInterpolationSamples;}
                }

                // Génère les triangles avec vertices dédupliqués.
                // Ordre 0, 2, 1 (pas 0, 1, 2) pour le winding attendu par RealtimeMesh.
                const uint64 StreamBuildingStart = bSampleCell
                    ? FPlatformTime::Cycles64() : 0;
                for (int32 i = 0; MCTriTable[CaseIndex][i] != -1; i += 3)
                {
                    const int32 E0 = MCTriTable[CaseIndex][i];
                    const int32 E1 = MCTriTable[CaseIndex][i + 1];
                    const int32 E2 = MCTriTable[CaseIndex][i + 2];
                    if (HybridResult.UsesDCEdge(cx, cy, cz, E0, CellsPerAxis))
                    {
                        continue;
                    }
                    const int32 Idx0 = GetOrCreateVertex(
                        EdgeVertices[E0], EdgeNormals[E0], EdgeSplitKeys[E0]);
                    const int32 Idx1 = GetOrCreateVertex(
                        EdgeVertices[E1], EdgeNormals[E1], EdgeSplitKeys[E1]);
                    const int32 Idx2 = GetOrCreateVertex(
                        EdgeVertices[E2], EdgeNormals[E2], EdgeSplitKeys[E2]);

                    TArray<int32>* Dst = &GroundTris;
                    if (!bStandaloneSkipMaterialGroups)
                    {
                        // F17 — majority vote (at least two ceiling vertices => ceiling triangle).
                        const int32 CapVotes = (int32)VertexClasses[Idx0]
                                             + (int32)VertexClasses[Idx1]
                                             + (int32)VertexClasses[Idx2];
                        Dst = (CapVotes >= 2) ? &CapTris : &GroundTris;
                    }
                    Dst->Add(Idx0);
                    Dst->Add(Idx2);
                    Dst->Add(Idx1);
                }
                if (bSampleCell)
                {
                    const uint64 StageCycles=FPlatformTime::Cycles64()-StreamBuildingStart;
                    if(bProfileMesher)VoxelDensityProfile::AddSampledMeasurement(
                        VoxelDensityProfile::EBucket::MesherStreamBuilding,StageCycles);
                    if(bAuditMesherStages)
                    {AuditDedupStreamCycles+=StageCycles;++AuditDedupStreamSamples;}
                }
            }
        }
    }

    // Add DC geometry to the same material runs as MC. Analytic feature IDs split corner normals.
    // Ajouter la geometrie DC dans les memes groupes que MC; les IDs separent les normales nettes.
    for (const VoxelHybridSharpMesher::FTriangle& Triangle : HybridResult.Triangles)
    {
        const int32 Idx0 = GetOrCreateVertex(
            Triangle.Positions[0], Triangle.Normals[0], Triangle.SplitKeys[0]);
        const int32 Idx1 = GetOrCreateVertex(
            Triangle.Positions[1], Triangle.Normals[1], Triangle.SplitKeys[1]);
        const int32 Idx2 = GetOrCreateVertex(
            Triangle.Positions[2], Triangle.Normals[2], Triangle.SplitKeys[2]);
        const int32 CapVotes = bStandaloneSkipMaterialGroups ? 0
            : static_cast<int32>(VertexClasses[Idx0])
                + static_cast<int32>(VertexClasses[Idx1])
                + static_cast<int32>(VertexClasses[Idx2]);
        TArray<int32>& Dst = (CapVotes >= 2) ? CapTris : GroundTris;
        Dst.Add(Idx0);
        Dst.Add(Idx2);
        Dst.Add(Idx1);
    }

    if (bProfileMesher)
    {
        VoxelDensityProfile::AddMeasurement(
            VoxelDensityProfile::EBucket::MesherCellClassification,
            0, CellClassificationCalls);
        VoxelDensityProfile::AddMeasurement(
            VoxelDensityProfile::EBucket::MesherGradientNormals,
            0, SurfaceCellCalls);
        VoxelDensityProfile::AddMeasurement(
            VoxelDensityProfile::EBucket::MesherVertexInterpolation,
            0, SurfaceCellCalls);
        VoxelDensityProfile::AddMeasurement(
            VoxelDensityProfile::EBucket::MesherStreamBuilding,
            0, SurfaceCellCalls);
    }

    VoxelDensityProfile::FScopedTimer MesherFinalisationTimer(
        VoxelDensityProfile::EBucket::MesherOther);

    //=========================================================================
    // SKIRTS — boucheurs de fissures aux coutures de LOD (clipmap)
    //=========================================================================
    // Deux tuiles de niveaux voisins maillent à des résolutions différentes : leurs iso-surfaces
    // ne se rejoignent pas parfaitement le long de la face partagée → une fine fissure traversante.
    // On la scelle en suspendant une courte « jupe » (mur) sous chaque arête de surface posée sur
    // l'une des 6 faces externes de la tuile, extrudée VERS LE SOLIDE le long de la normale inversée
    // sur ~une taille de cellule. Là où la surface du voisin est décalée, les jupes des deux tuiles
    // se recouvrent dans la roche et ferment le trou ; ailleurs la jupe est enterrée et invisible.
    // Émis en DOUBLE FACE (deux orientations) pour s'afficher quel que soit le côté caméra / le
    // matériau. Une arête de surface posée sur une face de frontière a sa coordonnée d'axe EXACTE
    // (l'interpolation MC garde fixe l'axe de la face) → comparaison flottante exacte fiable.
    MeshData.NumSurfaceVertices = MeshData.Vertices.Num();
    const int32 GroundSurfaceTriangleCount=GroundTris.Num()/3;
    const int32 CapSurfaceTriangleCount=CapTris.Num()/3;
    TArray<uint8> GroundSkirtFaceMasks,CapSkirtFaceMasks;
    const uint64 SkirtStageStartCycles=OutStageTimings?FPlatformTime::Cycles64():0;
    l77MarchingCubesStage.Restore();
    VoxelForgeCore::Linux77Profiler::FScopedStage l77SkirtStage(
        VoxelForgeCore::Linux77Profiler::EStage::Skirt);
    VoxelForgeCore::floating_island_detail::FScopedIslandFieldCaller skirtFieldCaller(
        VoxelForgeCore::floating_island_detail::EIslandFieldCaller::Skirt);
    if (bGenerateSkirts && (GroundTris.Num() + CapTris.Num()) > 0)
    {
        const float ExtentCm = (float)(CellsPerAxis * Step) * GeometryVoxelSize;
        const float MinX = ChunkWorldPos.X, MinY = ChunkWorldPos.Y, MinZ = ChunkWorldPos.Z;
        const float MaxX = MinX + ExtentCm, MaxY = MinY + ExtentCm, MaxZ = MinZ + ExtentCm;
        const float SkirtDepth = FMath::Max(1.0f, SkirtCells) * (float)Step * GeometryVoxelSize;

        auto BoundaryFaceMask = [&](const FVector& A, const FVector& B) -> uint8
        {
            uint8 Mask=0;
            if(A.X==MinX&&B.X==MinX)Mask|=1u<<0;
            if(A.X==MaxX&&B.X==MaxX)Mask|=1u<<1;
            if(A.Y==MinY&&B.Y==MinY)Mask|=1u<<2;
            if(A.Y==MaxY&&B.Y==MaxY)Mask|=1u<<3;
            if(A.Z==MinZ&&B.Z==MinZ)Mask|=1u<<4;
            if(A.Z==MaxZ&&B.Z==MaxZ)Mask|=1u<<5;
            return Mask;
        };

        auto AddSkirtVert = [&](const FVector& Pos, const FVector& Nrm, const FColor& Col) -> int32
        {
            const int32 Idx = MeshData.Vertices.Num();
            MeshData.Vertices.Add(Pos);
            MeshData.Normals.Add(Nrm);
            if (bEmitUV)
            {
                MeshData.UVs.Add(FVector2f(Pos.X / VOXEL_SIZE, Pos.Y / VOXEL_SIZE));
            }
            if (bEmitColors)
            {
                MeshData.Colors.Add(Col);   // parallel with the other vertex arrays
            }
            return Idx;
        };

        std::vector<float> SurfaceVertexDensitySamples;
        const auto CountSkirtDensityProbe=[&](EVoxelSkirtDensityProbeCaller Caller,
            double X,double Y,double Z)
        {
            if(!OutStageTimings)return;
            const std::size_t CallerIndex=static_cast<std::size_t>(Caller);
            ++OutStageTimings->SkirtDensityProbes;
            ++OutStageTimings->SkirtDensityProbesByCaller[CallerIndex];
            const auto IsLatticeCoordinate=[](double Value)
            {
                const double QueryValue=static_cast<double>(static_cast<float>(Value));
                return std::isfinite(QueryValue)
                    &&std::abs(QueryValue-std::round(QueryValue))<=1.0e-5;
            };
            if(!IsLatticeCoordinate(X)||!IsLatticeCoordinate(Y)
                ||!IsLatticeCoordinate(Z))
                ++OutStageTimings->SkirtDensityOffGridProbesByCaller[CallerIndex];
        };
        const auto QuerySkirtDensity=[&](EVoxelSkirtDensityProbeCaller Caller,
            double X,double Y,double Z)
        {
            CountSkirtDensityProbe(Caller,X,Y,Z);
            EPassageDensityQueryCaller QueryCaller=
                EPassageDensityQueryCaller::SkirtSurfaceVertex;
            switch(Caller)
            {
                case EVoxelSkirtDensityProbeCaller::SurfaceVertex:
                    QueryCaller=EPassageDensityQueryCaller::SkirtSurfaceVertex;break;
                case EVoxelSkirtDensityProbeCaller::ExtrudedVertex:
                    QueryCaller=EPassageDensityQueryCaller::SkirtExtrudedVertex;break;
                case EVoxelSkirtDensityProbeCaller::ClippedVertex:
                    QueryCaller=EPassageDensityQueryCaller::SkirtClippedVertex;break;
                case EVoxelSkirtDensityProbeCaller::RefinementFace:
                    QueryCaller=EPassageDensityQueryCaller::SkirtRefinementFace;break;
                case EVoxelSkirtDensityProbeCaller::CoarseCentroid:
                    QueryCaller=EPassageDensityQueryCaller::SkirtCoarseCentroid;break;
                case EVoxelSkirtDensityProbeCaller::Count:break;
            }
            FStandaloneShapeWorld::FScopedPassageDensityQueryCaller queryCaller(QueryCaller);
            return ShapeWorldContext->DensityAt(X,Y,Z);
        };
        if(ShapeWorldContext)
            SurfaceVertexDensitySamples.assign(static_cast<size_t>(MeshData.Vertices.Num()),
                std::numeric_limits<float>::quiet_NaN());
        const auto SampleSurfaceVertexDensity=[&](int32 OriginalIndex,const FVector& Position)
        {
            if(!ShapeWorldContext)return 0.0f;
            if(OriginalIndex<0||static_cast<size_t>(OriginalIndex)>=SurfaceVertexDensitySamples.size())
            {
                const double X=Position.X/GeometryVoxelSize;
                const double Y=Position.Y/GeometryVoxelSize;
                const double Z=Position.Z/GeometryVoxelSize;
                return QuerySkirtDensity(EVoxelSkirtDensityProbeCaller::SurfaceVertex,
                    X,Y,Z);
            }
            float& Density=SurfaceVertexDensitySamples[static_cast<size_t>(OriginalIndex)];
            if(std::isnan(Density))
            {
                const double X=Position.X/GeometryVoxelSize;
                const double Y=Position.Y/GeometryVoxelSize;
                const double Z=Position.Z/GeometryVoxelSize;
                Density=QuerySkirtDensity(EVoxelSkirtDensityProbeCaller::SurfaceVertex,
                    X,Y,Z);
            }
            return Density;
        };

        // F17 — jupes émises PAR SEAU : chaque jupe hérite la classe (donc le polygroup /
        // matériau) de son triangle source ; ajouter les jupes après coup casserait les runs
        // d'indices contigus par groupe qu'exige RMC.
        auto EmitSkirts = [&](TArray<int32>& Tris,TArray<uint8>& FaceMasks)
        {
            // On ajoute en itérant : on fige le nombre de triangles de surface et on n'ajoute qu'au-delà.
            const int32 BaseTriNum = Tris.Num();
            struct FSkirtVertex
            {
                FVector Position;
                FVector Normal;
                FColor Color;
                float Density=0.0f;
                int32 OriginalIndex=INDEX_NONE;
            };
            const auto MakeSkirtVertex=[&](const FVector& Position,const FVector& Normal,
                const FColor& Color,bool bOnSurface,int32 OriginalIndex)
            {
                FSkirtVertex Vertex{Position,Normal,Color,0.0f,OriginalIndex};
                if(!bOnSurface&&ShapeWorldContext)
                {
                    const double X=Position.X/GeometryVoxelSize;
                    const double Y=Position.Y/GeometryVoxelSize;
                    const double Z=Position.Z/GeometryVoxelSize;
                    Vertex.Density=QuerySkirtDensity(
                        EVoxelSkirtDensityProbeCaller::ExtrudedVertex,X,Y,Z);
                }
                return Vertex;
            };
            const auto InterpolateVertex=[&](const FSkirtVertex& A,const FSkirtVertex& B,float T)
            {
                if(T<=1.0e-6f)return A;
                if(T>=1.0f-1.0e-6f)return B;
                FSkirtVertex Out;
                Out.Position=FMath::Lerp(A.Position,B.Position,T);
                Out.Normal=FMath::Lerp(A.Normal,B.Normal,T).GetSafeNormal();
                Out.Color=FColor(
                    static_cast<uint8>(FMath::Clamp(FMath::RoundToInt(FMath::Lerp(static_cast<float>(A.Color.R),static_cast<float>(B.Color.R),T)),0,255)),
                    static_cast<uint8>(FMath::Clamp(FMath::RoundToInt(FMath::Lerp(static_cast<float>(A.Color.G),static_cast<float>(B.Color.G),T)),0,255)),
                    static_cast<uint8>(FMath::Clamp(FMath::RoundToInt(FMath::Lerp(static_cast<float>(A.Color.B),static_cast<float>(B.Color.B),T)),0,255)),
                    static_cast<uint8>(FMath::Clamp(FMath::RoundToInt(FMath::Lerp(static_cast<float>(A.Color.A),static_cast<float>(B.Color.A),T)),0,255)));
                if(ShapeWorldContext)
                {
                    const double X=Out.Position.X/GeometryVoxelSize;
                    const double Y=Out.Position.Y/GeometryVoxelSize;
                    const double Z=Out.Position.Z/GeometryVoxelSize;
                    Out.Density=QuerySkirtDensity(
                        EVoxelSkirtDensityProbeCaller::ClippedVertex,X,Y,Z);
                }
                else Out.Density=FMath::Lerp(A.Density,B.Density,T);
                Out.OriginalIndex=INDEX_NONE;
                return Out;
            };
            const auto EmitClippedTriangle=[&](auto&& Self,const FSkirtVertex& A,
                const FSkirtVertex& B,const FSkirtVertex& C,uint8 FaceMask,
                int32 RefinementDepth)->void
            {
                if(GVoxelForgeRawSkirtsForAudit)
                {
                    const FSkirtVertex Input[3]={A,B,C};
                    int32 Indices[3];
                    for(int32 Vertex=0;Vertex<3;++Vertex)
                        Indices[Vertex]=Input[Vertex].OriginalIndex>=0
                            ?Input[Vertex].OriginalIndex
                            :AddSkirtVert(Input[Vertex].Position,Input[Vertex].Normal,
                                Input[Vertex].Color);
                    Tris.Add(Indices[0]);Tris.Add(Indices[1]);Tris.Add(Indices[2]);
                    Tris.Add(Indices[0]);Tris.Add(Indices[2]);Tris.Add(Indices[1]);
                    FaceMasks.Add(FaceMask);FaceMasks.Add(FaceMask);
                    return;
                }
                FSkirtVertex Input[3]={A,B,C};
                FSkirtVertex Polygon[5];
                int32 PolygonCount=0;
                for(int32 Edge=0;Edge<3;++Edge)
                {
                    const FSkirtVertex& Start=Input[Edge];
                    const FSkirtVertex& End=Input[(Edge+1)%3];
                    const bool bStartSolid=Start.Density<=0.0f;
                    const bool bEndSolid=End.Density<=0.0f;
                    if(bStartSolid!=bEndSolid)
                    {
                        const float Denominator=Start.Density-End.Density;
                        const float T=FMath::Abs(Denominator)>1.0e-8f
                            ?FMath::Clamp(Start.Density/Denominator,0.0f,1.0f):0.5f;
                        Polygon[PolygonCount++]=InterpolateVertex(Start,End,T);
                    }
                    if(bEndSolid)Polygon[PolygonCount++]=End;
                }
                if(PolygonCount<3)return;
                int32 PolygonIndices[5];
                for(int32 Vertex=0;Vertex<PolygonCount;++Vertex)
                    PolygonIndices[Vertex]=Polygon[Vertex].OriginalIndex>=0
                        ?Polygon[Vertex].OriginalIndex
                        :AddSkirtVert(Polygon[Vertex].Position,Polygon[Vertex].Normal,
                            Polygon[Vertex].Color);
                for(int32 Fan=1;Fan+1<PolygonCount;++Fan)
                {
                    const FSkirtVertex& P0=Polygon[0];
                    const FSkirtVertex& P1=Polygon[Fan];
                    const FSkirtVertex& P2=Polygon[Fan+1];
                    if(ShapeWorldContext&&Step<=2)
                    {
                        const auto ProbeDensity=[&](double WA,double WB,double WC)
                        {
                            const FVector Position=P0.Position*static_cast<float>(WA)
                                +P1.Position*static_cast<float>(WB)
                                +P2.Position*static_cast<float>(WC);
                            const double X=Position.X/GeometryVoxelSize;
                            const double Y=Position.Y/GeometryVoxelSize;
                            const double Z=Position.Z/GeometryVoxelSize;
                            return QuerySkirtDensity(
                                EVoxelSkirtDensityProbeCaller::RefinementFace,X,Y,Z);
                        };
                        const auto VertexProbeDensity=[&](const FSkirtVertex& Vertex)
                        {
                            return Vertex.OriginalIndex>=0
                                ?SampleSurfaceVertexDensity(Vertex.OriginalIndex,Vertex.Position)
                                :Vertex.Density;
                        };
                        float MaximumProbeDensity=FMath::Max(VertexProbeDensity(P0),
                            FMath::Max(VertexProbeDensity(P1),VertexProbeDensity(P2)));
                        const double ProbeWeights[4][3]={
                            {1.0/3.0,1.0/3.0,1.0/3.0},
                            {0.5,0.5,0.0},{0.0,0.5,0.5},{0.5,0.0,0.5}};
                        const int32 ProbeCount=RefinementDepth==0?4:1;
                        for(int32 Probe=0;Probe<ProbeCount;++Probe)
                        {
                            MaximumProbeDensity=FMath::Max(MaximumProbeDensity,
                                ProbeDensity(ProbeWeights[Probe][0],
                                    ProbeWeights[Probe][1],ProbeWeights[Probe][2]));
                            if(MaximumProbeDensity>0.25f)break;
                        }
                        if(MaximumProbeDensity>0.25f)
                        {
                            if(RefinementDepth>=1)
                            {
                                continue;
                            }
                            const auto MakeMidpoint=[&](const FSkirtVertex& X,
                                const FSkirtVertex& Y)
                            {
                                FSkirtVertex Mid=InterpolateVertex(X,Y,0.5f);
                                Mid.OriginalIndex=INDEX_NONE;
                                return Mid;
                            };
                            const FSkirtVertex M01=MakeMidpoint(P0,P1);
                            const FSkirtVertex M12=MakeMidpoint(P1,P2);
                            const FSkirtVertex M20=MakeMidpoint(P2,P0);
                            Self(Self,P0,M01,M20,FaceMask,RefinementDepth+1);
                            Self(Self,M01,P1,M12,FaceMask,RefinementDepth+1);
                            Self(Self,M20,M12,P2,FaceMask,RefinementDepth+1);
                            Self(Self,M01,M12,M20,FaceMask,RefinementDepth+1);
                            continue;
                        }
                    }
                    else if(ShapeWorldContext)
                    {
                        const FVector Center=(P0.Position+P1.Position+P2.Position)
                            /(3.0f*GeometryVoxelSize);
                        if(QuerySkirtDensity(EVoxelSkirtDensityProbeCaller::CoarseCentroid,
                            Center.X,Center.Y,Center.Z)>0.0f)continue;
                    }
                    Tris.Add(PolygonIndices[0]);Tris.Add(PolygonIndices[Fan]);Tris.Add(PolygonIndices[Fan+1]);
                    Tris.Add(PolygonIndices[0]);Tris.Add(PolygonIndices[Fan+1]);Tris.Add(PolygonIndices[Fan]);
                    FaceMasks.Add(FaceMask);FaceMasks.Add(FaceMask);
                }
            };
            for (int32 t = 0; t + 2 < BaseTriNum; t += 3)
            {
                const int32 Tri[3] = { Tris[t], Tris[t + 1], Tris[t + 2] };
                for (int32 e = 0; e < 3; ++e)
                {
                    if(OutStageTimings)++OutStageTimings->SkirtEdgeCandidates;
                    const int32 iA = Tri[e], iB = Tri[(e + 1) % 3];
                    // COPIES par valeur — AddSkirtVert réalloue Vertices/Normals (invaliderait des refs).
                    const FVector PA = MeshData.Vertices[iA];
                    const FVector PB = MeshData.Vertices[iB];
                    const uint8 FaceMask=BoundaryFaceMask(PA,PB);
                    if (!FaceMask) continue;
                    if(OutStageTimings)++OutStageTimings->SkirtBoundaryCandidates;

                    const FVector NA = MeshData.Normals[iA];
                    const FVector NB = MeshData.Normals[iB];
                    const FColor  CA = bEmitColors ? MeshData.Colors[iA] : FColor::White;
                    const FColor  CB = bEmitColors ? MeshData.Colors[iB] : FColor::White;
                    const FSkirtVertex A=MakeSkirtVertex(PA,NA,CA,true,iA);
                    const FSkirtVertex B=MakeSkirtVertex(PB,NB,CB,true,iB);
                    const FSkirtVertex A2=MakeSkirtVertex(PA-NA*SkirtDepth,NA,CA,false,INDEX_NONE);
                    const FSkirtVertex B2=MakeSkirtVertex(PB-NB*SkirtDepth,NB,CB,false,INDEX_NONE);

                    // Clip both source triangles against the solid half-space before emitting
                    // each orientation. A curtain that points into tunnel air must not become
                    // visible geometry beyond the sampled isosurface.
                    EmitClippedTriangle(EmitClippedTriangle,A,B,B2,FaceMask,0);
                    EmitClippedTriangle(EmitClippedTriangle,A,B2,A2,FaceMask,0);
                }
            }
            MeshData.NumSkirtTriangles += (Tris.Num() - BaseTriNum) / 3;
        };
            EmitSkirts(GroundTris,GroundSkirtFaceMasks);
            EmitSkirts(CapTris,CapSkirtFaceMasks);
        }
    if(OutStageTimings)
        OutStageTimings->SkirtMilliseconds=static_cast<double>(
            FPlatformTime::Cycles64()-SkirtStageStartCycles)
                *FPlatformTime::GetSecondsPerCycle64()*1000.0;
    l77SkirtStage.Restore();
    VoxelForgeCore::Linux77Profiler::FScopedStage l77MesherFinalStage(
        VoxelForgeCore::Linux77Profiler::EStage::MarchingCubes);

    // F17 — concatène sol PUIS sky-cap : un run d'indices contigu par polygroup (voir
    // FVoxelMeshData::NumCeilingTriangles). Mêmes triangles qu'avant, seul l'ordre change.
    MeshData.Triangles.Reserve(GroundTris.Num() + CapTris.Num());
    MeshData.Triangles.Append(GroundTris);
    MeshData.Triangles.Append(CapTris);
    MeshData.NumCeilingTriangles = CapTris.Num() / 3;
    MeshData.TriangleSkirtFaceMasks.Reset();
    MeshData.TriangleSkirtFaceMasks.Reserve(MeshData.Triangles.Num()/3);
    for(int32 Triangle=0;Triangle<GroundSurfaceTriangleCount;++Triangle)
        MeshData.TriangleSkirtFaceMasks.Add(0);
    MeshData.TriangleSkirtFaceMasks.Append(GroundSkirtFaceMasks);
    for(int32 Triangle=0;Triangle<CapSurfaceTriangleCount;++Triangle)
        MeshData.TriangleSkirtFaceMasks.Add(0);
    MeshData.TriangleSkirtFaceMasks.Append(CapSkirtFaceMasks);

    if (VoxelDensityProfile::AreCountersEnabled())
    {
        const uint64 OutputArrayBytes = static_cast<uint64>(MeshData.Vertices.Num()) * sizeof(FVector)
            + static_cast<uint64>(MeshData.Normals.Num()) * sizeof(FVector)
            + static_cast<uint64>(MeshData.UVs.Num()) * sizeof(FVector2D)
            + static_cast<uint64>(MeshData.Colors.Num()) * sizeof(FColor)
            + static_cast<uint64>(MeshData.Triangles.Num()) * sizeof(int32);
        const uint64 OutputAllocatedBytes = static_cast<uint64>(MeshData.Vertices.GetAllocatedSize())
            + static_cast<uint64>(MeshData.Normals.GetAllocatedSize())
            + static_cast<uint64>(MeshData.UVs.GetAllocatedSize())
            + static_cast<uint64>(MeshData.Colors.GetAllocatedSize())
            + static_cast<uint64>(MeshData.Triangles.GetAllocatedSize());
        VoxelDensityProfile::AddCounter(
            VoxelDensityProfile::ECounter::MesherOutputArrayBytes,
            OutputArrayBytes);
        VoxelDensityProfile::AddCounter(
            VoxelDensityProfile::ECounter::MesherOutputAllocatedBytes,
            OutputAllocatedBytes);
    }

    if(OutStageTimings)
    {
        OutStageTimings->CellClassificationCalls=CellClassificationCalls;
        OutStageTimings->SurfaceCellCalls=SurfaceCellCalls;
        const double CycleToMs=FPlatformTime::GetSecondsPerCycle64()*1000.0;
        const auto ScaleSampledCycles=[CycleToMs](uint64 Cycles,uint64 Samples,uint64 Calls)
        {
            return Samples>0?static_cast<double>(Cycles)*static_cast<double>(Calls)
                /static_cast<double>(Samples)*CycleToMs:0.0;
        };
        OutStageTimings->CellClassificationMilliseconds=
            static_cast<double>(AuditCellClassificationCycles)*CycleToMs;
        OutStageTimings->GradientNormalsMilliseconds=ScaleSampledCycles(
            AuditGradientCycles,AuditGradientSamples,SurfaceCellCalls);
        OutStageTimings->EdgeInterpolationMilliseconds=ScaleSampledCycles(
            AuditInterpolationCycles,AuditInterpolationSamples,SurfaceCellCalls);
        OutStageTimings->VertexDedupStreamMilliseconds=ScaleSampledCycles(
            AuditDedupStreamCycles,AuditDedupStreamSamples,SurfaceCellCalls);
        OutStageTimings->SkirtMilliseconds=std::max(0.0,
            OutStageTimings->SkirtMilliseconds);
        const double TotalMilliseconds=static_cast<double>(
            FPlatformTime::Cycles64()-StageMesherStartCycles)*CycleToMs;
        OutStageTimings->OtherMilliseconds=FMath::Max(0.0,TotalMilliseconds
            -OutStageTimings->DensityGridMilliseconds
            -OutStageTimings->CellClassificationMilliseconds
            -OutStageTimings->GradientNormalsMilliseconds
            -OutStageTimings->EdgeInterpolationMilliseconds
            -OutStageTimings->VertexDedupStreamMilliseconds
            -OutStageTimings->SkirtMilliseconds);
    }

    return MeshData;
}