Mesh Segmentation

I need to segment a mesh for laser processing. In areas with high Gaussian curvature, I plan to divide the mesh into small regions, while in areas with low Gaussian curvature, I want to create larger blocks to reduce machine movement and improve laser processing efficiency. Currently, I’m using the k-means method to segment the mesh, which can segment it into blocks of specified dimensions, but I don’t know how to merge the regions with low Gaussian curvature. Could someone help me? Additionally, I have both surface and mesh data available, which can both be used @dale
MeshSegmentation.3dm (964.3 KB)

this is what I want


this is what I get

this is my code

public class MeshSegmentation
{
    /// <summary>
    /// </summary>
    public static void SegmentMeshBySize(
        Mesh inputMesh,
        double maxSize,
        out List<Mesh> segmentedMeshes,
        out List<Color> colors)
    {
        segmentedMeshes = new List<Mesh>();
        colors = new List<Color>();

        var faceIndices = new List<int>();
        for (int i = 0; i < inputMesh.Faces.Count; i++)
            faceIndices.Add(i);

        var results = RecursiveSegment(faceIndices, inputMesh, maxSize, 2);

        segmentedMeshes.AddRange(results);

        var rnd = new Random();
        foreach (var m in segmentedMeshes)
        {
            Color color = Color.FromArgb(255, rnd.Next(256), rnd.Next(256), rnd.Next(256));
            colors.Add(color);
        }
    }

    /// <summary>
    /// </summary>
    private static List<Mesh> RecursiveSegment(
        List<int> faceIndices,
        Mesh sourceMesh,
        double maxSize,
        int clusterCount)
    {
        var outputMeshes = new List<Mesh>();
        if (faceIndices.Count == 0)
            return outputMeshes;

        var faceCenters = new List<Point3d>();
        foreach (int i in faceIndices)
            faceCenters.Add(GetFaceCenter(sourceMesh, sourceMesh.Faces[i]));

        int[] labels = KMeansCluster(faceCenters, clusterCount);

        var clusters = new Dictionary<int, List<int>>();
        for (int i = 0; i < labels.Length; i++)
        {
            int label = labels[i];
            if (!clusters.ContainsKey(label))
                clusters[label] = new List<int>();
            clusters[label].Add(faceIndices[i]);
        }

        foreach (var kvp in clusters)
        {
            var indices = kvp.Value;
            Mesh m = ExtractMeshFromFaces(sourceMesh, indices);
            if (m == null || !m.IsValid)
                continue;

            var bbox = m.GetBoundingBox(true);
            double dx = bbox.Max.X - bbox.Min.X;
            double dy = bbox.Max.Y - bbox.Min.Y;
            double dz = bbox.Max.Z - bbox.Min.Z;

            if (dx > maxSize || dy > maxSize || dz > maxSize)
            {
                var subMeshes = RecursiveSegment(indices, sourceMesh, maxSize, 2);
                outputMeshes.AddRange(subMeshes);
            }
            else
            {
                outputMeshes.Add(m);
            }
        }

        return outputMeshes;
    }

    /// <summary>
    /// </summary>
    private static Mesh ExtractMeshFromFaces(Mesh source, List<int> faceIndices)
    {
        Mesh result = new Mesh();
        var vertMap = new Dictionary<int, int>();

        foreach (int i in faceIndices)
        {
            var face = source.Faces[i];
            int[] indices = face.IsQuad
                ? new[] { face.A, face.B, face.C, face.D }
                : new[] { face.A, face.B, face.C };

            foreach (int vi in indices)
            {
                if (!vertMap.ContainsKey(vi))
                {
                    int newIndex = result.Vertices.Count;
                    result.Vertices.Add(source.Vertices[vi]);
                    vertMap[vi] = newIndex;
                }
            }

            if (face.IsQuad)
                result.Faces.AddFace(vertMap[face.A], vertMap[face.B], vertMap[face.C], vertMap[face.D]);
            else
                result.Faces.AddFace(vertMap[face.A], vertMap[face.B], vertMap[face.C]);
        }

        result.Normals.ComputeNormals();
        result.Compact();
        return result;
    }

    /// <summary>
    /// </summary>
    private static Point3d GetFaceCenter(Mesh mesh, MeshFace face)
    {
        var a = mesh.Vertices[face.A];
        var b = mesh.Vertices[face.B];
        var c = mesh.Vertices[face.C];

        if (face.IsQuad)
        {
            var d = mesh.Vertices[face.D];
            return new Point3d((a.X + b.X + c.X + d.X) / 4.0,
                               (a.Y + b.Y + c.Y + d.Y) / 4.0,
                               (a.Z + b.Z + c.Z + d.Z) / 4.0);
        }
        else
        {
            return new Point3d((a.X + b.X + c.X) / 3.0,
                               (a.Y + b.Y + c.Y) / 3.0,
                               (a.Z + b.Z + c.Z) / 3.0);
        }
    }

    /// <summary>
    /// </summary>
    private static int[] KMeansCluster(List<Point3d> points, int k, int maxIterations = 100)
    {
        var rand = new Random();
        int n = points.Count;
        var centroids = new List<Point3d>();
        var labels = new int[n];

        for (int i = 0; i < k; i++)
            centroids.Add(points[rand.Next(n)]);

        for (int iter = 0; iter < maxIterations; iter++)
        {
            for (int i = 0; i < n; i++)
            {
                double minDist = double.MaxValue;
                int best = 0;
                for (int j = 0; j < k; j++)
                {
                    double dist = points[i].DistanceTo(centroids[j]);
                    if (dist < minDist)
                    {
                        minDist = dist;
                        best = j;
                    }
                }
                labels[i] = best;
            }

            var newCentroids = new Point3d[k];
            var counts = new int[k];

            for (int j = 0; j < k; j++)
                newCentroids[j] = new Point3d(0, 0, 0);

            for (int i = 0; i < n; i++)
            {
                int label = labels[i];
                newCentroids[label] += points[i];
                counts[label]++;
            }

            for (int j = 0; j < k; j++)
                if (counts[j] > 0)
                    newCentroids[j] = newCentroids[j] / counts[j];

            centroids = new List<Point3d>(newCentroids);
        }

        return labels;
    }
}