pmneila / PyMCubes

Marching cubes (and related tools) for Python
BSD 3-Clause "New" or "Revised" License
692 stars 87 forks source link

Interpolation for Vertex #43

Closed denghilbert closed 1 year ago

denghilbert commented 1 year ago

I noticed this lib implemented the interpolation for vertex like this:

double mc_isovalue_interpolation(double isovalue, double f1, double f2, double x1, double x2)
{
    if(f2==f1)
        return (x2+x1)/2;

    return (x2-x1)*(isovalue-f1)/(f2-f1) + x1;
}

But when I run:

field = np.ones([2, 2, 2]) + 99
field[0][0][0] = 0

vertices1, triangles1 = mcubes.marching_cubes(field, 50)

init = np.ones([2, 2, 2]) + 99
init[0][0][0] = 40

vertices2, triangles2 = mcubes.marching_cubes(field, 50)

The result is like:

vertices1
array([[0.5, 0. , 0. ],
       [0. , 0.5, 0. ],
       [0. , 0. , 0.5]])
vertices2
array([[0.5, 0. , 0. ],
       [0. , 0.5, 0. ],
       [0. , 0. , 0.5]])

It seems that it does not interpolate based on the value of the corner. If there is any bug or I just misunderstood this. :sob:

pmneila commented 1 year ago

Hi @denghilbert,

There is a typo in your code 😅. You are defining two different arrays, field and init, but then you use only field when calling marching_cubes. The fixed code works as expected:

field = np.ones([2, 2, 2]) + 99
field[0][0][0] = 0

vertices1, triangles1 = mcubes.marching_cubes(field, 50)

# Change `init` to `field`
field = np.ones([2, 2, 2]) + 99
field[0][0][0] = 40

vertices2, triangles2 = mcubes.marching_cubes(field, 50)

Output:

In [11]: vertices1
Out[11]:
array([[0.5, 0. , 0. ],
       [0. , 0.5, 0. ],
       [0. , 0. , 0.5]])

In [12]: vertices2
Out[12]:
array([[0.16666667, 0.        , 0.        ],
       [0.        , 0.16666667, 0.        ],
       [0.        , 0.        , 0.16666667]])

Hope that helps.

Best!