diff --git a/src/gpytoolbox/reach_for_the_arcs.py b/src/gpytoolbox/reach_for_the_arcs.py index 96a8ef11..f788c1ae 100644 --- a/src/gpytoolbox/reach_for_the_arcs.py +++ b/src/gpytoolbox/reach_for_the_arcs.py @@ -11,7 +11,7 @@ def reach_for_the_arcs(U, S, rng_seed=3452, return_point_cloud=False, - fine_tune_iters=10, + fine_tune_iters=3, batch_size=10000, num_rasterization_spheres=0, screening_weight=10., @@ -22,6 +22,7 @@ def reach_for_the_arcs(U, S, local_search_t=0.01, tol=1e-4, clamp_value=np.Inf, + force_cpu=False, parallel=False, verbose=False): """Creates a mesh from a signed distance function (SDF) using the @@ -69,6 +70,8 @@ def reach_for_the_arcs(U, S, tolerance for determining whether a point is inside a sphere clamp_value : float, optional (default np.Inf) value to which the SDF is clamped for clamped SDF reconstruction + force_cpu : bool, optional (default False) + whether to force rasterization onto the CPU parallel : bool, optional (default False) whether to parallelize the algorithm or not verbose : bool, optional (default False) @@ -137,6 +140,7 @@ def reach_for_the_arcs(U, S, local_search_iters=local_search_iters, batch_size=batch_size, num_rasterization_spheres=num_rasterization_spheres, + force_cpu=force_cpu, tol=tol, clamp_value=clamp_value, parallel=parallel, verbose=verbose) else: @@ -151,6 +155,7 @@ def reach_for_the_arcs(U, S, local_search_iters=local_search_iters, batch_size=batch_size, num_rasterization_spheres=num_rasterization_spheres, + force_cpu=force_cpu, tol=tol, clamp_value=clamp_value, parallel=parallel, verbose=verbose) else: @@ -171,6 +176,7 @@ def reach_for_the_arcs(U, S, n_local_searches=n_local_searches, local_search_iters=local_search_iters, batch_size=batch_size, + force_cpu=force_cpu, tol=tol, clamp_value=clamp_value, parallel=parallel, verbose=verbose) @@ -254,6 +260,7 @@ def _sdf_to_point_cloud(U, S, num_rasterization_spheres=0, tol=1e-4, clamp_value=np.Inf, + force_cpu=False, parallel=False, verbose=False): """Converts an SDF to a point cloud where all points are valid with respect @@ -327,7 +334,7 @@ def _sdf_to_point_cloud(U, S, P = _outside_points_from_rasterization(U, S, rng_seed=seed(), res=rasterization_resolution, num_spheres=num_rasterization_spheres, - tol=tol, parallel=parallel, verbose=verbose) + tol=tol, force_cpu=force_cpu, parallel=parallel, verbose=verbose) if verbose: print(f" Rasterization took {time.time()-t0_rasterization}s") @@ -664,10 +671,7 @@ def _fine_tune_point_cloud(U, S, P, N, f, tol, clamp_value, parallel, verbose) if verbose: - if debug_Vgt is not None and debug_Fgt is not None: - print(f" After fine tuning iter {it}, we have {f.size} points. (Current chamfer error = {chamfer(V, F, debug_Vgt, debug_Fgt)})") - else: - print(f" After fine tuning iter {it}, we have {f.size} points.") + print(f" After fine tuning iter {it}, we have {f.size} points.") return P, N, f diff --git a/test/test_reach_for_the_arcs.py b/test/test_reach_for_the_arcs.py index 97eb1039..bbfc33fb 100644 --- a/test/test_reach_for_the_arcs.py +++ b/test/test_reach_for_the_arcs.py @@ -6,35 +6,36 @@ class TestReachForTheArcs(unittest.TestCase): def test_beat_marching_cubes_low_res(self): meshes = ["R.npy", "bunny_oded.obj", "armadillo.obj"] - ns = [10, 20, 30] for mesh in meshes: - for n in ns: - if mesh[-3:]=="obj": - v, f = gpy.read_mesh("test/unit_tests_data/" + mesh) - elif mesh[-3:]=="npy": - data = np.load("test/unit_tests_data/" + mesh, - allow_pickle=True) - v = data[()]['V'] - f = data[()]['F'] - v = gpy.normalize_points(v) - - sdf = lambda x: gpy.signed_distance(x, v, f)[0] - - if mesh[-3:]=="obj": - gx, gy, gz = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) - GV = np.vstack((gx.flatten(), gy.flatten(), gz.flatten())).T - V_mc, F_mc = gpy.marching_cubes(sdf(GV), GV, n+1, n+1, n+1) - elif mesh[-3:]=="npy": - gx, gy = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) - GV = np.vstack((gx.flatten(), gy.flatten())).T - V_mc, F_mc = gpy.marching_squares(sdf(GV), GV, n+1, n+1) - - h_mc = gpy.approximate_hausdorff_distance(V_mc, F_mc.astype(np.int32), v, f.astype(np.int32), use_cpp = True) - U,G = gpy.reach_for_the_arcs(GV, sdf(GV), parallel=False) - h_ours = gpy.approximate_hausdorff_distance(U, G.astype(np.int32), v, f.astype(np.int32), use_cpp = True) - - # print(f"reach_for_the_arcs h: {h_ours}, MC h: {h_mc} for {mesh} with n={n}") - self.assertTrue(h_ours < h_mc) + if mesh[-3:]=="obj": + v, f = gpy.read_mesh("test/unit_tests_data/" + mesh) + elif mesh[-3:]=="npy": + data = np.load("test/unit_tests_data/" + mesh, + allow_pickle=True) + v = data[()]['V'] + f = data[()]['F'] + v = gpy.normalize_points(v) + + sdf = lambda x: gpy.signed_distance(x, v, f)[0] + n = 10 + + if mesh[-3:]=="obj": + gx, gy, gz = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) + GV = np.vstack((gx.flatten(), gy.flatten(), gz.flatten())).T + V_mc, F_mc = gpy.marching_cubes(sdf(GV), GV, n+1, n+1, n+1) + elif mesh[-3:]=="npy": + gx, gy = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) + GV = np.vstack((gx.flatten(), gy.flatten())).T + V_mc, F_mc = gpy.marching_squares(sdf(GV), GV, n+1, n+1) + + h_mc = gpy.approximate_hausdorff_distance(V_mc, F_mc.astype(np.int32), v, f.astype(np.int32), use_cpp = True) + U,G = gpy.reach_for_the_arcs(GV, sdf(GV), fine_tune_iters=3, + local_search_iters=3, + parallel=True, verbose=False) + h_ours = gpy.approximate_hausdorff_distance(U, G.astype(np.int32), v, f.astype(np.int32), use_cpp = True) + + #print(f"reach_for_the_arcs h: {h_ours}, MC h: {h_mc} for {mesh} with n={n}") + self.assertTrue(h_ours < h_mc) def test_noop(self): @@ -50,7 +51,7 @@ def test_noop(self): v = gpy.normalize_points(v) sdf = lambda x: gpy.signed_distance(x, v, f)[0] - n = 20 + n = 10 if mesh[-3:]=="obj": gx, gy, gz = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) @@ -59,35 +60,65 @@ def test_noop(self): gx, gy = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) GV = np.vstack((gx.flatten(), gy.flatten())).T - U,G = gpy.reach_for_the_arcs(GV, sdf(GV), parallel=False) - + U,G = gpy.reach_for_the_arcs(GV, sdf(GV), fine_tune_iters=3, + local_search_iters=3, + parallel=True, verbose=False) + h = gpy.approximate_hausdorff_distance(U, G.astype(np.int32), v, f.astype(np.int32), use_cpp=True) self.assertTrue(h < 0.1) - Up,Gp = gpy.reach_for_the_arcs(GV, sdf(GV), parallel=True) + + def test_parallel_is_the_same(self): + meshes = ["R.npy", "bunny_oded.obj", "armadillo.obj"] + for mesh in meshes: + if mesh[-3:]=="obj": + v, f = gpy.read_mesh("test/unit_tests_data/" + mesh) + elif mesh[-3:]=="npy": + data = np.load("test/unit_tests_data/" + mesh, + allow_pickle=True) + v = data[()]['V'] + f = data[()]['F'] + v = gpy.normalize_points(v) + + sdf = lambda x: gpy.signed_distance(x, v, f)[0] + n = 10 + + if mesh[-3:]=="obj": + gx, gy, gz = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) + GV = np.vstack((gx.flatten(), gy.flatten(), gz.flatten())).T + elif mesh[-3:]=="npy": + gx, gy = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) + GV = np.vstack((gx.flatten(), gy.flatten())).T + + U,G = gpy.reach_for_the_arcs(GV, sdf(GV), fine_tune_iters=3, + local_search_iters=3, + parallel=False, verbose=False) + + Up,Gp = gpy.reach_for_the_arcs(GV, sdf(GV), fine_tune_iters=3, + local_search_iters=3, + parallel=True, verbose=False) h_parallel = gpy.approximate_hausdorff_distance(U, G.astype(np.int32), Up, Gp.astype(np.int32), use_cpp = True) # print(f"parallel Hausdorff distance h: {h_parallel}") self.assertTrue(h_parallel < 1e-6) + def test_simple_is_sdf_violated(self): - meshes = ["cube.obj", "hemisphere.obj", "cone.obj"] + meshes = ["cube.obj", "hemisphere.obj"] for mesh in meshes: - n = 30 + n = 10 v, f = gpy.read_mesh("test/unit_tests_data/" + mesh) v = gpy.normalize_points(v) sdf = lambda x: gpy.signed_distance(x, v, f)[0] gx, gy, gz = np.meshgrid(np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1), np.linspace(-1.0, 1.0, n+1)) GV = np.vstack((gx.flatten(), gy.flatten(), gz.flatten())).T - U,G = gpy.reach_for_the_arcs(GV, sdf(GV), parallel=False) + U,G = gpy.reach_for_the_arcs(GV, sdf(GV), fine_tune_iters=3, + local_search_iters=3, + parallel=True, verbose=False) sdf_rec = lambda x: gpy.signed_distance(x, U, G)[0] - self.assertTrue(np.max(np.abs(sdf(GV)-sdf_rec(GV))) < 0.02) - - Up,Gp = gpy.reach_for_the_arcs(GV, sdf(GV), parallel=True) - h_parallel = gpy.approximate_hausdorff_distance(U, G.astype(np.int32), Up, Gp.astype(np.int32), use_cpp = True) - # print(f"parallel Hausdorff distance h: {h_parallel}") - self.assertTrue(h_parallel < 1e-6) + print(np.max(np.abs(sdf(GV)-sdf_rec(GV)))) + self.assertTrue(np.max(np.abs(sdf(GV)-sdf_rec(GV))) < 0.03) if __name__ == '__main__':