52

I'm working on a Marching Cubes implementation in Unity. My code is based on Paul Bourke's code actually with a lot of modifications, but anyway I'm checking if a block at a position is null if it is than a debug texture will be placed on it.

Here is an image of the problem

This is my MC script.

    public class MarchingCubes
{
private World world;
private Chunk chunk;
private List<Vector3> vertices = new List<Vector3> ();
private List<Vector3> normals = new List<Vector3> ();
private Vector3[] ns;
private List<int> triangles = new List<int> ();
private List<Vector2> uvs = new List<Vector2> ();
private Vector3[] positions = new Vector3[8];
private float[] corners = new float[8];
private Vector3i size = new Vector3i (16, 128, 16);

Vector3[] vertlist = new Vector3[12];

private float isolevel = 1f;

private float Corner (Vector3i pos)
{
    int x = pos.x;
    int y = pos.y;
    int z = pos.z;
    if (x < size.x && z < size.z) {
        return chunk.GetValue (x, y, z);
    } else {
        int ix = chunk.X, iz = chunk.Z;
        int rx = chunk.region.x, rz = chunk.region.z;
        if (x >= size.x) {
            ix++;
            x = 0;
        }

        if (z >= size.z) {
            iz++;
            z = 0;
        }
        return chunk.region.GetChunk (ix, iz).GetValue (x, y, z);
    }
}

Block block;

public Mesh MarchChunk (World world, Chunk chunk, Mesh mesh)
{
    this.world = world;
    this.chunk = chunk;

    vertices.Clear ();
    triangles.Clear ();
    uvs.Clear ();

    for (int x = 0; x < size.x; x++) {
        for (int y = 1; y < size.y - 2; y++) {
            for (int z = 0; z < size.z; z++) {

                block = chunk.GetBlock (x, y, z);
                int cubeIndex = 0;

                for (int i = 0; i < corners.Length; i++) {
                    corners [i] = Corner (new Vector3i (x, y, z) + offset [i]);
                    positions [i] = (new Vector3i (x, y, z) + offset [i]).ToVector3 ();

                    if (corners [i] < isolevel)
                        cubeIndex |= (1 << i);
                }

                if (eTable [cubeIndex] == 0)
                    continue;

                for (int i = 0; i < vertlist.Length; i++) {
                    if ((eTable [cubeIndex] & 1 << i) == 1 << i)
                        vertlist [i] = LinearInt (positions [eCons [i, 0]], positions [eCons [i, 1]], corners [eCons [i, 0]], corners [eCons [i, 1]]);
                }

                for (int i = 0; triTable [cubeIndex, i] != -1; i += 3) {
                    int index = vertices.Count;

                    vertices.Add (vertlist [triTable [cubeIndex, i]]);
                    vertices.Add (vertlist [triTable [cubeIndex, i + 1]]);
                    vertices.Add (vertlist [triTable [cubeIndex, i + 2]]);

                    float tec = (0.125f);
                    Vector2 uvBase = block != null ? block.UV : new Vector2 ();

                    uvs.Add (uvBase);
                    uvs.Add (uvBase + new Vector2 (0, tec));
                    uvs.Add (uvBase + new Vector2 (tec, tec));
                                    
                    triangles.Add (index + 0);
                    triangles.Add (index + 1);
                    triangles.Add (index + 2);
                }
            }
        }
    }       

    if (mesh == null)
        mesh = new Mesh ();
    mesh.Clear ();
    mesh.vertices = vertices.ToArray ();
    mesh.triangles = triangles.ToArray ();
    mesh.uv = uvs.ToArray ();
    mesh.RecalculateNormals ();
    return mesh;
}

bool IsBitSet (int b, int pos)
{
    return ((b & pos) == pos);
}

Vector3 LinearInt (Vector3 p1, Vector3 p2, float v1, float v2)
{
    Vector3 p;
    p.x = p1.x + (isolevel - v1) * (p2.x - p1.x) / (v2 - v1);
    p.y = p1.y + (isolevel - v1) * (p2.y - p1.y) / (v2 - v1);
    p.z = p1.z + (isolevel - v1) * (p2.z - p1.z) / (v2 - v1);
    return p;
}

private static int[,] eCons = new int[12, 2] {
    { 0, 1 },
    { 1, 2 },
    { 2, 3 },
    { 3, 0 },
    { 4, 5 },
    { 5, 6 },
    { 6, 7 },
    { 7, 4 },
    { 0, 4 },
    { 1, 5 },
    { 2, 6 },
    { 3, 7 }
};

private static Vector3i[] offset = new Vector3i[8] {
    new Vector3i (0, 0, 1),
    new Vector3i (1, 0, 1),
    new Vector3i (1, 0, 0),
    new Vector3i (0, 0, 0),
    new Vector3i (0, 1, 1),
    new Vector3i (1, 1, 1),
    new Vector3i (1, 1, 0),
    new Vector3i (0, 1, 0)
};
  }

I didn't put the tables in the sample, because they are the same as the ones in Bourke's code.

EDIT: What I figured out yet is that the cell's value at the blue triangles is 0 so they don't have to be triangulated, but the cell's value under them is 1 and because of this a top triangle is created to complete the mesh.

HarrY
  • 607
  • 4
  • 14
Statey
  • 631
  • 5
  • 10
  • 6
    It looks a table problem, because it's always the same `cubeIndex`. I have seen implementations that bit 6 and 7 are swapped in this table: `if (grid.val[6] < isolevel) cubeindex |= 128; if (grid.val[7] < isolevel) cubeindex |= 64;` Look here at this line: [WEBGL_Marchingcubes of threejs.org](https://github.com/mrdoob/three.js/blob/master/examples/js/MarchingCubes.js#L179) _(I've noticed this because I'm working on a implementation on the GPU using a geometryshader)_ – Jeroen van Langen Jun 26 '17 at 12:59
  • 8
    Perhaps not relevant but here goes: The original algorithm has a flaw that can result in holes in the generated mesh. Many years ago, I wrote [a paper about this](https://martinliversage.blob.core.windows.net/publications/1992%20Automedica.pdf). I also wrote a master’s thesis (in Danish) where I touched upon the issue. I’ve tried to illustrate in on [page 81 in my thesis](https://martinliversage.blob.core.windows.net/publications/1994-09%20Visualisering%20af%20tredimensionale%20medicinske%20data.pdf). On page 87 you can see all possible meshes where the troublesome are marked with an *. – Martin Liversage Jun 26 '17 at 13:15
  • 2
    Okay I'll try this and as soon as I get something I'll tell you about it. – Statey Jun 26 '17 at 13:20
  • So @JeroenvanLangen this is interresting but it gives me errors because of the non-finite values at specific vertex indexes – Statey Jun 26 '17 at 13:27
  • @MartinLiversage I'll have a look at it – Statey Jun 26 '17 at 13:27
  • What does the `chunk.GetValue(...)` returns? – Jeroen van Langen Jun 26 '17 at 13:49
  • Ohh it's a shortcut it's the same as `chunk.GetBlock(x,y,z).GetValue()` – Statey Jun 26 '17 at 13:58
  • It gets the value of the given cell. – Statey Jun 26 '17 at 14:22
  • @JeroenvanLangen do you have any ideas? – Statey Jun 26 '17 at 14:49
  • @MartinLiversage I had a look at you thesis and the troublesome meshes, but how can I fixed the problems with the given cases? – Statey Jun 26 '17 at 15:20
  • 4
    @Statey: I'm not sure if the problems that your experience are a result of the ambiguities in the original algorithm so I can't give you specific advice. However, a strategy for handling the problem of cube meshes with holes is to perform interpolation of the center of the ambiguous side of a cube (or the center of an ambiguous cube) to determine if the interpolated point should be inside or outside the isosurface. This is then used to create the correct isosurface mesh and because the decision is made deterministically for all cubes the mesh will be connected without holes. – Martin Liversage Jun 27 '17 at 13:35
  • In the `Corner` function, try to check for negative as well `if (x >= 0 && x < size.x && z >= 0 && z < size.z)`, then handle accordingly in the `else` block – Sorashi Dec 08 '22 at 16:35
  • Requires more input from OP to solve, and OP's last answer was 6 years ago. The problem isn't general enough that we can hope for an answer that'd be useful to anyone else. – hugo Aug 31 '23 at 07:37

0 Answers0