From 5cbb0d31ac307c4aeb38761c7aa4c4b60db54bef Mon Sep 17 00:00:00 2001 From: gdamms Date: Wed, 1 Feb 2023 09:12:18 +0100 Subject: [PATCH] good enough --- src/draw.py | 51 +++++++++++++------- src/sfs2d.py | 134 ++++++++++++++++++++++++++++++++++----------------- 2 files changed, 124 insertions(+), 61 deletions(-) diff --git a/src/draw.py b/src/draw.py index 36ec1b0..261ba38 100644 --- a/src/draw.py +++ b/src/draw.py @@ -7,8 +7,8 @@ VOXEL_SIZE = 0.1 X_MIN, X_MAX = -5, 5 Y_MIN, Y_MAX = -5, 5 -x_vals = np.arange(X_MIN, X_MAX, VOXEL_SIZE) -y_vals = np.arange(Y_MIN, Y_MAX, VOXEL_SIZE) +x_vals = np.arange(X_MIN + VOXEL_SIZE / 2, X_MAX, VOXEL_SIZE) +y_vals = np.arange(Y_MIN + VOXEL_SIZE / 2, Y_MAX, VOXEL_SIZE) image_length = 1500 @@ -23,15 +23,15 @@ Z = f(X - 2, Y) + f(X + 2, Y) Z = (Z > 0.4).astype(np.float32) Z *= np.random.rand(*Z.shape) -for i, x in enumerate(x_vals): - for j, y in enumerate(y_vals): - color = f"{hex(int(Z[j, i] * 255))[2:]}" +for j, x in enumerate(x_vals): + for i, y in enumerate(y_vals): + color = f"{hex(int(Z[i, j] * 255))[2:]}" if color == "0": color = "#003" else: color = "#" + 3 * color - plt.fill([x, x + 0.1, x + 0.1, x], - [y, y, y + 0.1, y + 0.1], color=color) + plt.fill([x - VOXEL_SIZE/2, x + VOXEL_SIZE/2, x + VOXEL_SIZE/2, x - VOXEL_SIZE/2], + [y - VOXEL_SIZE/2, y - VOXEL_SIZE/2, y + VOXEL_SIZE/2, y + VOXEL_SIZE/2], color=color) nb_cams = 32 cam_poses = np.array( @@ -49,8 +49,8 @@ cam2world_projs = np.array( ) for i in range(nb_cams): - x = np.array([[0, 0, 1], [0.5, -0.2, 1], [0.5, 0.2, 1], [0, 0, 1]]).T - x = cam2world_projs[i] @ x + x = np.array([[0, 0, 1], [0.5, -0.2, 1], [0.5, 0.2, 1], [0, 0, 1]]) + x = cam2world_projs[i] @ x.T plt.plot(cam_poses[i][0], cam_poses[i][1], "ro") plt.text(cam_poses[i][0], cam_poses[i][1], str(i)) plt.plot(x[0, :], x[1, :], "r-") @@ -79,16 +79,15 @@ for i in range(nb_cams): for j in pixels_sort: x, y = X.flatten()[j], Y.flatten()[j] - color = f"{hex(int(Z.flatten()[j] * 255))[2:]}" - if color == "0": + color = int(Z.flatten()[j] * 255) + if color == 0: continue - color = "#" + 3 * color RT = np.linalg.inv(cam2world_projs[i])[:-1, :] - px = np.array([[x, y, 1], - [x + 0.1, y, 1], - [x + 0.1, y + 0.1, 1], - [x, y + 0.1, 1]]) + px = np.array([[x - VOXEL_SIZE/2, y - VOXEL_SIZE/2, 1], + [x + VOXEL_SIZE/2, y - VOXEL_SIZE/2, 1], + [x + VOXEL_SIZE/2, y + VOXEL_SIZE/2, 1], + [x - VOXEL_SIZE/2, y + VOXEL_SIZE/2, 1]]) px = RT @ px.T px /= px[0, :] px += np.array([[0], [1.0]]) @@ -98,8 +97,26 @@ for i in range(nb_cams): x1 = px[1, :].max() img[:, int(x0):int(x1), :] = np.array( - 3 * [int(Z.flatten()[j] * 255)], dtype=np.uint8) + 3 * [color], dtype=np.uint8) mask[:, int(x0):int(x1), :] = np.array([255, 255, 255], dtype=np.uint8) + # x, y = X_MIN + 27 * VOXEL_SIZE, Y_MIN + 66 * VOXEL_SIZE + # RT = np.linalg.inv(cam2world_projs[i])[:-1, :] + # px = np.array([[x + VOXEL_SIZE/2, y + VOXEL_SIZE/2, 1], + # [x, y, 1], + # [x + VOXEL_SIZE, y, 1], + # [x + VOXEL_SIZE, y + VOXEL_SIZE, 1], + # [x, y + VOXEL_SIZE, 1]]) + # px = RT @ px.T + # px /= px[0, :] + # px += np.array([[0], [1.0]]) + # px *= 0.5 * np.array([[1], [image_length]]) + + # img[48:52, int(px[1, 0]-1):int(px[1, 0]+1), :] = np.array([255, 0, 0], dtype=np.uint8) + # img[48:52, int(px[1, 1]-1):int(px[1, 1]+1), :] = np.array([255, 255, 0], dtype=np.uint8) + # img[48:52, int(px[1, 2]-1):int(px[1, 2]+1), :] = np.array([0, 255, 255], dtype=np.uint8) + # img[48:52, int(px[1, 3]-1):int(px[1, 3]+1), :] = np.array([255, 0, 255], dtype=np.uint8) + # img[48:52, int(px[1, 4]-1):int(px[1, 4]+1), :] = np.array([0, 255, 0], dtype=np.uint8) + cv2.imwrite(f"data/peanut/images/Image{i:04}.png", img) cv2.imwrite(f"data/peanut/masks/Image{i:04}.png", mask) diff --git a/src/sfs2d.py b/src/sfs2d.py index 9e1b083..a954d2f 100644 --- a/src/sfs2d.py +++ b/src/sfs2d.py @@ -8,14 +8,11 @@ from borders import update_border from fvi import fast_voxel_intersect - -x_vals = np.linspace(-5, 5, 100) -y_vals = np.linspace(-5, 5, 100) - VOXEL_SIZE = 0.1 X_MIN, X_MAX = -5, 5 Y_MIN, Y_MAX = -5, 5 + def plot_camera(cam2world_proj, name, color): points = np.array([[0, 0, 1], [0.5, -0.2, 1], [0.5, 0.2, 1], [0, 0, 1]]).T points = cam2world_proj @ points @@ -26,12 +23,12 @@ def plot_camera(cam2world_proj, name, color): def plot_voxels(voxels, color, alpha, background): if background: - plt.fill([X_MIN, X_MAX, X_MAX, X_MIN], [Y_MIN, Y_MIN, Y_MAX, Y_MAX], color="#000") + plt.fill([X_MIN, X_MAX, X_MAX, X_MIN], [ + Y_MIN, Y_MIN, Y_MAX, Y_MAX], color="#000") for i, j in np.argwhere(voxels): x, y = X_MIN + i * VOXEL_SIZE, Y_MIN + j * VOXEL_SIZE - plt.fill([x, x + VOXEL_SIZE, x + VOXEL_SIZE, x], [y, y, y + VOXEL_SIZE, y + VOXEL_SIZE], color=color, alpha=alpha) - - + plt.fill([x, x + VOXEL_SIZE, x + VOXEL_SIZE, x], [y, y, y + + VOXEL_SIZE, y + VOXEL_SIZE], color=color, alpha=alpha) nb_frame = 32 @@ -40,8 +37,8 @@ points = np.array( [ [x, y, 1.0] for x, y in product( - np.arange(X_MIN, X_MAX, VOXEL_SIZE), - np.arange(Y_MIN, Y_MAX, VOXEL_SIZE), + np.arange(X_MIN + VOXEL_SIZE / 2, X_MAX, VOXEL_SIZE), + np.arange(Y_MIN + VOXEL_SIZE / 2, Y_MAX, VOXEL_SIZE), ) ] ) @@ -54,7 +51,8 @@ cam_rots = np.linspace(np.pi, 3 * np.pi, nb_frame, endpoint=False) cam2world_projs = np.array( [ [[np.cos(theta), -np.sin(theta), cam_pose[0]], - [np.sin(theta), np.cos(theta), cam_pose[1]], [0, 0, 1]] + [np.sin(theta), np.cos(theta), cam_pose[1]], + [0, 0, 1]] for theta, cam_pose in zip(cam_rots, cam_poses) ] ) @@ -68,18 +66,16 @@ frames = [] for k in range(nb_frame): frame = cv2.imread( f"data/peanut/masks/Image{k:04}.png", cv2.IMREAD_GRAYSCALE)[0, :] - frames.append(cv2.imread(f"data/peanut/images/Image{k:04}.png", cv2.IMREAD_GRAYSCALE)[0, :]) - # print(f"frame.shape: {frame.shape}") + frames.append(cv2.imread( + f"data/peanut/images/Image{k:04}.png", cv2.IMREAD_GRAYSCALE)[0, :]) RT = np.linalg.inv(cam2world_projs[k])[:-1, :] - # a = np.array([0.0, 0.0, 1.0]) cam_points = RT @ points.T cam_points /= cam_points[0] cam_points += np.array([[0], [1.0]]) cam_points *= 0.5 * np.array([[1], [frame.shape[0]]]) cam_points = np.round(cam_points).astype(np.int32) - # print(f"cam_points: {cam_points}") visible = np.logical_and.reduce( ( @@ -94,8 +90,10 @@ for k in range(nb_frame): cam_points = cam_points[:, solid] points = points[solid, :] -voxel = np.zeros((int((X_MAX-X_MIN)/VOXEL_SIZE), int((Y_MAX-Y_MIN)/VOXEL_SIZE))) -idx = np.floor_divide(points[:, :2] - np.array([X_MIN, Y_MIN]), VOXEL_SIZE).astype(int) +voxel = np.zeros((int((X_MAX-X_MIN)/VOXEL_SIZE), + int((Y_MAX-Y_MIN)/VOXEL_SIZE))) +idx = np.floor_divide( + points[:, :2] - np.array([X_MIN, Y_MIN]), VOXEL_SIZE).astype(int) voxel[idx[:, 0], idx[:, 1]] = 1 @@ -103,15 +101,14 @@ def f(x, y): return np.exp(-((x**2) + y**2) / 3) -x_vals = np.linspace(-5, 5, 100) -y_vals = np.linspace(-5, 5, 100) - +x_vals = np.arange(X_MIN + VOXEL_SIZE / 2, X_MAX, VOXEL_SIZE) +y_vals = np.arange(Y_MIN + VOXEL_SIZE / 2, Y_MAX, VOXEL_SIZE) X, Y = np.meshgrid(x_vals, y_vals) Z = f(X - 2, Y) + f(X + 2, Y) Z = (Z > 0.4).astype(np.float32) -Z = np.swapaxes(Z, 0, 1) +Z = Z.T plt.figure() plt.axis("equal") @@ -119,37 +116,60 @@ plot_voxels(Z, "#fff", 1.0, True) plot_voxels(voxel, "#00f", 0.5, False) for i in range(nb_frame): plot_camera(cam2world_projs[i], str(i), "#f00") - +# plt.show() +plt.savefig(f"init.png", dpi=300, bbox_inches="tight", pad_inches=0, transparent=True) +plt.close() origin = np.array([X_MIN, Y_MIN]) step = np.array([VOXEL_SIZE, VOXEL_SIZE]) -shape = np.array([int((X_MAX-X_MIN)/VOXEL_SIZE), int((Y_MAX-Y_MIN)/VOXEL_SIZE)]) +shape = np.array([int((X_MAX-X_MIN)/VOXEL_SIZE), + int((Y_MAX-Y_MIN)/VOXEL_SIZE)]) +# x, y = 2.45, -1.55 +# i, j = np.floor_divide(np.array([x, y]) - origin, step).astype(int) +# print(i, j) +# exit() -for _ in range(5): +for p in range(11): border = update_border(voxel) + voxel_ = voxel.copy() plt.figure() plt.axis("equal") - plot_voxels(voxel, "#fff", 1.0, True) + plot_voxels(Z, "#fff", 1.0, True) # plot_voxels(border, "#00f", 0.5, False) for i in range(nb_frame): plot_camera(cam2world_projs[i], str(i), "#f00") + # plt.show() for i, j in track(np.argwhere(border)): - x, y = X_MIN + i * VOXEL_SIZE, Y_MIN + j * VOXEL_SIZE - start = np.array([x, y]) + 0.5 * step - - # plt.fill([x, x + VOXEL_SIZE, x + VOXEL_SIZE, x], [y, y, y + VOXEL_SIZE, y + VOXEL_SIZE], "#0f0", alpha=0.7) + # if i != 74 or j != 34: + # continue + x, y = X_MIN + (i + 0.5) * VOXEL_SIZE, Y_MIN + (j + 0.5) * VOXEL_SIZE + start = np.array([x, y]) + # plt.plot([start[0]], [start[1]], "o", color="#f00") + # plt.fill([x, x + VOXEL_SIZE, x + VOXEL_SIZE, x], + # [y, y, y + VOXEL_SIZE, y + VOXEL_SIZE], "#0f0", alpha=0.7) values = [] + # a b c d all four corners + a = np.array([x, y]) + # plt.plot([a[0]], [a[1]], "o", color="#ff0") + b = np.array([x + VOXEL_SIZE, y]) + # plt.plot([b[0]], [b[1]], "o", color="#0ff") + c = np.array([x + VOXEL_SIZE, y + VOXEL_SIZE]) + # plt.plot([c[0]], [c[1]], "o", color="#f0f") + d = np.array([x, y + VOXEL_SIZE]) + # plt.plot([d[0]], [d[1]], "o", color="#0f0") + for k in range(nb_frame): # plot_camera(cam2world_projs[k], str(k), "#0f0") end = cam_poses[k] - _, _, voxels_intersected = fast_voxel_intersect(start, end, origin, step, shape) + _, _, voxels_intersected = fast_voxel_intersect( + start, end, origin, step, shape) voxels_intersected = np.array(voxels_intersected, dtype=np.int32) # for voxel_id in voxels_intersected: @@ -175,29 +195,49 @@ for _ in range(5): if len(voxels_intersected) == 0: visible = True else: - visible = voxel[voxels_intersected[:, 0], voxels_intersected[:, 1]].sum() == 0 + visible = voxel[voxels_intersected[:, 0], + voxels_intersected[:, 1]].sum() == 0 if visible: RT = np.linalg.inv(cam2world_projs[k])[:-1, :] - proj = RT @ np.array([start[0], start[1], 1.0]) - proj /= proj[0] - proj += np.array([0, 1.0]) - proj *= 0.5 * np.array([1, frames[k].shape[0]]) + # concatenate start with a b c and d + points = np.concatenate([start[np.newaxis, :], + a[np.newaxis, :], + b[np.newaxis, :], + c[np.newaxis, :], + d[np.newaxis, :]], axis=0) + # homogenous coordinates + points = np.concatenate([points, np.ones((5, 1))], axis=1) + proj = RT @ points.T + proj /= proj[0, :] + proj += np.array([[0], [1.0]]) + proj *= 0.5 * np.array([[1], [frames[k].shape[0]]]) proj = np.round(proj).astype(np.int32) - values.append(frames[k][proj[1]]) + values.append(frames[k][proj[1, 0]]) - # plt.figure() - # plt.imshow(np.repeat(frames[k][np.newaxis, :], 100, axis=0), cmap="gray") - # plt.plot(proj[1], 50, "ro") - # plt.show() + # fig = plt.gcf() + # plt.figure() + # plt.imshow( + # np.repeat(frames[k][np.newaxis, :], 100, axis=0), cmap="gray") + # plt.plot(proj[1, 0], 50, "o", color="#f00") + # plt.plot(proj[1, 1], 50, "o", color="#ff0") + # plt.plot(proj[1, 2], 50, "o", color="#0ff") + # plt.plot(proj[1, 3], 50, "o", color="#f0f") + # plt.plot(proj[1, 4], 50, "o", color="#0f0") + # plt.figure(fig) # else: # plot_camera(cam2world_projs[k], str(k), "#f00") - if np.std(values) < 2.0: - voxel[i, j] = 1 + # print(values) + # plt.show() + + occurences = [[x, values.count(x)] for x in set(values)] + occurences = sorted(occurences, key=lambda x: x[1], reverse=True) + if occurences[0][1] >= len(values) * 0.68: + voxel_[i, j] = 1 color = "#0f0" else: - voxel[i, j] = 0 + voxel_[i, j] = 0 color = "#f00" plt.fill( [ @@ -216,9 +256,15 @@ for _ in range(5): alpha=0.3, ) + voxel = voxel_ + + plt.savefig(f"evol{p}.png", dpi=300, bbox_inches="tight", pad_inches=0, transparent=True) + plt.figure() plt.axis("equal") plot_voxels(Z, "#fff", 1.0, True) plot_voxels(voxel, "#00f", 0.5, False) + plt.savefig(f"shape{p}.png", dpi=300, bbox_inches="tight", pad_inches=0, transparent=True) + plt.close("all") - plt.show() + # plt.show()