Issue in Unsigned Distance Field for open Meshes

[pinakinathc]

Hi, I am trying to calculate the unsigned distance field of an open mesh using:

signed_distance, new_faces, new_verts = igl.signed_distance(samples, mesh.vertices, mesh.faces, return_normals=False)
unsigned_distance = signed_distance**2

But the generated mesh gives a wrong value near the open region -- it behaves as if the mesh was closed (i.e., watertight).

I am attaching the code to recreate the issue. Is this a problem of the package?

Sample data: https://drive.google.com/file/d/1rkfDW5cne-0x8AFeVvT-9tDkiBDZQNnL/view?usp=sharing

# -*- coding: utf-8 -*-
# author: pinakinathc

import os
import igl
import numpy as np
import trimesh
from skimage import measure


def create_grid(resX, resY, resZ, b_min=np.array([0, 0, 0]), b_max=np.array([1, 1, 1]), transform=None):
    '''
    Create a dense grid of given resolution and bounding box
    :param resX: resolution along X axis
    :param resY: resolution along Y axis
    :param resZ: resolution along Z axis
    :param b_min: vec3 (x_min, y_min, z_min) bounding box corner
    :param b_max: vec3 (x_max, y_max, z_max) bounding box corner
    :return: [3, resX, resY, resZ] coordinates of the grid, and transform matrix from mesh index
    '''
    coords = np.mgrid[:resX, :resY, :resZ]
    coords = coords.reshape(3, -1)
    coords_matrix = np.eye(4)
    length = b_max - b_min
    coords_matrix[0, 0] = length[0] / resX
    coords_matrix[1, 1] = length[1] / resY
    coords_matrix[2, 2] = length[2] / resZ
    coords_matrix[0:3, 3] = b_min
    coords = np.matmul(coords_matrix[:3, :3], coords) + coords_matrix[:3, 3:4]
    if transform is not None:
        coords = np.matmul(transform[:3, :3], coords) + transform[:3, 3:4]
        coords_matrix = np.matmul(transform, coords_matrix)
    coords = coords.reshape(3, resX, resY, resZ)
    return coords, coords_matrix


def save_obj_mesh(mesh_path, verts, faces):
    file = open(mesh_path, 'w')

    print ('shape of vertices: {}, faces: {}'.format(
        verts.shape, faces.shape))
    for v in verts:
        file.write('v %.4f %.4f %.4f\n' % (v[0], v[1], v[2]))
    for f in faces:
        f_plus = f + 1
        file.write('f %d %d %d\n' % (f_plus[0], f_plus[2], f_plus[1]))
    file.close()


def save_samples_truncted_prob(fname, points, prob):
    ''' Save the visualization of sampling to a ply file '''
    red = prob.reshape([-1, 1])/prob.max() * 255
    green = red
    blue = np.zeros(red.shape)

    to_save = np.concatenate([points, red, green, blue], axis=-1)
    return np.savetxt(fname, to_save, fmt='%.6f %.6f %.6f %d %d %d', comments='',
        header=('ply\nformat ascii 1.0\nelement vertex {:d}\nproperty float x\nproperty float y\
            \nproperty float z\nproperty uchar red\nproperty uchar green\nproperty uchar blue\nend_header').format(points.shape[0]))


if __name__ == '__main__':
	datapath = 'shirt_mesh_r_tmp.obj'
	mesh = trimesh.load(datapath)
	## Rescale
	scene_obj = trimesh.Scene(mesh)
	scene_obj = scene_obj.scaled(0.9/scene_obj.scale)
	mesh = scene_obj.geometry['shirt_mesh_r_tmp.obj']
	mesh.rezero()
	vertices = mesh.vertices
	b_max = vertices.max(0)
	b_min = vertices.min(0)
	delta = 0.2

	resolution = 64
	samples, calibs = create_grid(resolution, resolution, resolution, b_min-delta, b_max+delta)
	samples = samples.reshape(3, -1).T # Nx3

	# Calculate Unsigned distance using IGL
	signed_distance, new_faces, new_verts = igl.signed_distance(samples, mesh.vertices, mesh.faces, return_normals=False)
	unsigned_distance = signed_distance**2

	# Generate mesh
	labels = unsigned_distance
	verts, faces, normals, values = measure.marching_cubes_lewiner(labels.reshape(resolution, resolution, resolution), 1e-4)
	verts = np.matmul(calibs[:3, :3], verts.T) + calibs[:3, 3:4]
	save_path = os.path.join(os.getcwd(), 'checking')
	save_obj_mesh(save_path+'.obj', verts.T, faces)
	save_samples_truncted_prob(save_path+'.ply', samples, labels)