clean up

.gitignore

@@ -14,6 +14,6 @@ orisys.egg-info/
 /models/
 SDK/
 
-/patch_generation/exported_data/
+/tools/3d_patch_generation_blender/exported_data
 /docs/exported_data/
 /docs/.quarto/

tools/3d_patch_generation_blender/blender_export.py

@@ -0,0 +1,147 @@
+import bpy
+import numpy as np
+import json
+import os
+from mathutils import Matrix, Vector
+
+def export_camera_intrinsics(camera_obj, scene):
+    """Return camera intrinsic matrix K (3x3) and distortion (none)."""
+    cam = camera_obj.data
+    render = scene.render
+    res_x = render.resolution_x
+    res_y = render.resolution_y
+
+    # Sensor size and focal length
+    if cam.sensor_fit == 'HORIZONTAL':
+        sensor_width = cam.sensor_width
+        sensor_height = sensor_width / (res_x / res_y)
+    else:
+        sensor_height = cam.sensor_height
+        sensor_width = sensor_height * (res_x / res_y)
+    focal = cam.lens  # in mm
+
+    # Pixel sizes in mm
+    pixel_width = sensor_width / res_x
+    pixel_height = sensor_height / res_y
+
+    # Intrinsic matrix (OpenCV style: fx, fy in pixels, cx, cy principal point)
+    fx = focal / pixel_width
+    fy = focal / pixel_height
+    cx = res_x / 2.0
+    cy = res_y / 2.0
+
+    K = np.array([[fx, 0, cx],
+                  [0, fy, cy],
+                  [0,  0,  1]], dtype=np.float64)
+    return K, res_x, res_y
+
+def export_mesh_data(obj, filepath):
+    """Export mesh vertices, faces, and UVs to a compressed .npz file."""
+    me = obj.data
+    # Transform vertices to world space
+    vertices = np.array([obj.matrix_world @ v.co for v in me.vertices], dtype=np.float32)
+    faces = np.array([poly.vertices[:] for poly in me.polygons], dtype=np.int32)
+
+    # UVs (per loop, need mapping to vertices for simplicity we'll store per-face)
+    uv_layer = me.uv_layers.active
+    if uv_layer:
+        uv_data = np.array([uv_layer.data[i].uv for i in range(len(uv_layer.data))], dtype=np.float32)
+        loop_to_vertex = np.array([me.loops[i].vertex_index for i in range(len(me.loops))], dtype=np.int32)
+    else:
+        uv_data = np.zeros((0,2), dtype=np.float32)
+        loop_to_vertex = np.zeros(0, dtype=np.int32)
+
+    np.savez_compressed(filepath,
+                        vertices=vertices,
+                        faces=faces,
+                        uv_data=uv_data,
+                        loop_to_vertex=loop_to_vertex)
+    print(f"Mesh exported to {filepath}")
+
+def export_camera_pose(camera_obj, filepath):
+    """Export camera extrinsics (world-to-camera transform matrix)."""
+    # Camera world matrix: transforms from camera local to world
+    cam_to_world = np.array(camera_obj.matrix_world, dtype=np.float32)
+    # We need world-to-camera transform (inverse)
+    world_to_cam = np.linalg.inv(cam_to_world)
+    np.save(filepath, world_to_cam)
+    print(f"Camera pose exported to {filepath}")
+
+def export_pixel_to_3d_lookup(obj, camera_obj, scene, filepath):
+    """
+    Precompute for every pixel the 3D hit point (world coordinates)
+    and save as a 2D array (height, width, 3). This avoids raycasting later.
+    """
+    depsgraph = bpy.context.evaluated_depsgraph_get()
+    render = scene.render
+    res_x = render.resolution_x
+    res_y = render.resolution_y
+
+    # Get camera data
+    cam = camera_obj.data
+    if cam.sensor_fit == 'HORIZONTAL':
+        sensor_width = cam.sensor_width
+        sensor_height = sensor_width / (res_x / res_y)
+    else:
+        sensor_height = cam.sensor_height
+        sensor_width = sensor_height * (res_x / res_y)
+    focal = cam.lens
+    frame_width = sensor_width / focal
+    frame_height = sensor_height / focal
+
+    # Preallocate array (height, width, 3) for XYZ world coordinates
+    xyz_map = np.full((res_y, res_x, 3), np.nan, dtype=np.float32)
+
+    origin = camera_obj.matrix_world.translation
+    cam_rot = camera_obj.matrix_world.to_quaternion()
+
+    for y in range(res_y):
+        # Normalized sensor coordinate v (0=bottom, 1=top)
+        v_ndc = (y + 0.5) / res_y
+        sy = (v_ndc - 0.5) * frame_height
+        for x in range(res_x):
+            u_ndc = (x + 0.5) / res_x
+            sx = (u_ndc - 0.5) * frame_width
+
+            direction_cam = Vector((sx, sy, -1.0))
+            direction_world = cam_rot @ direction_cam
+            direction_world.normalize()
+
+            hit, location, normal, face_index, hit_obj, matrix = scene.ray_cast(
+                depsgraph, origin, direction_world)
+            if hit and hit_obj == obj:
+                xyz_map[y, x, :] = location
+
+    np.savez_compressed(filepath, xyz_map=xyz_map)
+    print(f"Pixel-to-3D lookup saved to {filepath} (shape {res_y}x{res_x})")
+
+def export_all():
+    # Configuration
+    output_dir = bpy.path.abspath("//exported_data/")
+    os.makedirs(output_dir, exist_ok=True)
+
+    obj = bpy.context.active_object   # or specify name
+    cam = bpy.data.objects.get("Camera")
+    if not cam:
+        cam = bpy.context.scene.camera
+    if not obj or not cam:
+        raise RuntimeError("Select mesh and ensure camera exists.")
+
+    # 1. Mesh data
+    export_mesh_data(obj, os.path.join(output_dir, "mesh.npz"))
+
+    # 2. Camera intrinsics (K matrix) and resolution
+    K, w, h = export_camera_intrinsics(cam, bpy.context.scene)
+    np.savez(os.path.join(output_dir, "camera_intrinsics.npz"), K=K, width=w, height=h)
+
+    # 3. Camera pose (world-to-camera)
+    export_camera_pose(cam, os.path.join(output_dir, "world_to_cam.npy"))
+
+    # 4. Optional: precomputed pixel-to-3D map (fastest)
+    export_pixel_to_3d_lookup(obj, cam, bpy.context.scene,
+                              os.path.join(output_dir, "pixel_to_xyz.npz"))
+
+    print("Export complete.")
+
+if __name__ == "__main__":
+    export_all()

tools/3d_patch_generation_blender/mapto3d.py

@@ -0,0 +1,177 @@
+import numpy as np
+import cv2
+import os
+
+def load_exported_data(data_dir):
+    """Load all exported files."""
+    # Mesh
+    mesh = np.load(os.path.join(data_dir, "mesh.npz"))
+    vertices = mesh['vertices']
+    faces = mesh['faces']
+    uv_data = mesh['uv_data']
+    loop_to_vertex = mesh['loop_to_vertex']
+
+    # Camera intrinsics
+    cam_data = np.load(os.path.join(data_dir, "camera_intrinsics.npz"))
+    K = cam_data['K']
+    width = int(cam_data['width'])
+    height = int(cam_data['height'])
+
+    # Camera pose (world-to-camera)
+    world_to_cam = np.load(os.path.join(data_dir, "world_to_cam.npy"))
+
+    # Precomputed pixel-to-3D map (optional)
+    xyz_map_path = os.path.join(data_dir, "pixel_to_xyz.npz")
+    if os.path.exists(xyz_map_path):
+        xyz_map = np.load(xyz_map_path)['xyz_map']
+    else:
+        xyz_map = None
+
+    return {
+        'vertices': vertices,
+        'faces': faces,
+        'uv_data': uv_data,
+        'loop_to_vertex': loop_to_vertex,
+        'K': K,
+        'width': width,
+        'height': height,
+        'world_to_cam': world_to_cam,
+        'xyz_map': xyz_map
+    }
+
+def ray_intersect_triangle(ray_origin, ray_dir, v0, v1, v2):
+    """Möller–Trumbore intersection algorithm. Returns t (distance) or None."""
+    epsilon = 1e-8
+    edge1 = v1 - v0
+    edge2 = v2 - v0
+    h = np.cross(ray_dir, edge2)
+    a = np.dot(edge1, h)
+    if abs(a) < epsilon:
+        return None
+    f = 1.0 / a
+    s = ray_origin - v0
+    u = f * np.dot(s, h)
+    if u < 0.0 or u > 1.0:
+        return None
+    q = np.cross(s, edge1)
+    v = f * np.dot(ray_dir, q)
+    if v < 0.0 or u + v > 1.0:
+        return None
+    t = f * np.dot(edge2, q)
+    if t > epsilon:
+        return t
+    return None
+
+def pixel_to_3d_ray(u, v, K, world_to_cam):
+    """
+    Given pixel coordinates (u,v) and camera intrinsics K,
+    return a ray in world coordinates: origin and direction.
+    """
+    # Invert K to get normalized image coordinates
+    fx = K[0,0]; fy = K[1,1]; cx = K[0,2]; cy = K[1,2]
+    x_norm = (u - cx) / fx
+    y_norm = (v - cy) / fy
+
+    # Ray in camera space: (x_norm, y_norm, 1.0)
+    ray_cam = np.array([x_norm, y_norm, 1.0], dtype=np.float64)
+
+    # Transform to world
+    cam_to_world = np.linalg.inv(world_to_cam)
+    R = cam_to_world[:3, :3]
+    t = cam_to_world[:3, 3]
+    ray_world_dir = R @ ray_cam
+    ray_world_dir /= np.linalg.norm(ray_world_dir)
+    ray_world_origin = t
+
+    return ray_world_origin, ray_world_dir
+
+def map_image_to_3d(image_path, data, use_precomputed=True):
+    """
+    Map each pixel of the image to a 3D world point.
+    If use_precomputed and xyz_map exists, return that directly.
+    Otherwise, perform raycasting on the mesh.
+    Returns:
+        points: (H, W, 3) array of 3D points (NaN for misses).
+        colors: (H, W, 3) array of RGB colors from image.
+    """
+    img = cv2.imread(image_path)
+    img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)  # to RGB
+    h, w = img.shape[:2]
+
+    if use_precomputed and data['xyz_map'] is not None:
+        # Use precomputed map (must match resolution)
+        if h != data['height'] or w != data['width']:
+            raise ValueError("Image resolution does not match precomputed map.")
+        points = data['xyz_map'].copy()
+        colors = img.astype(np.float32) / 255.0
+        return points, colors
+
+    # Otherwise, perform raycasting on the mesh
+    vertices = data['vertices']
+    faces = data['faces']
+    K = data['K']
+    world_to_cam = data['world_to_cam']
+
+    points = np.full((h, w, 3), np.nan, dtype=np.float32)
+    colors = img.astype(np.float32) / 255.0
+
+    # For speed, we can use a bounding volume hierarchy (BVH).
+    # Here we implement a simple brute-force over all faces (slow for large meshes).
+    # In practice, precomputing the xyz_map is strongly recommended.
+
+    for y in range(h):
+        for x in range(w):
+            origin, direction = pixel_to_3d_ray(x, y, K, world_to_cam)
+            best_t = np.inf
+            best_hit = None
+            for face in faces:
+                v0 = vertices[face[0]]
+                v1 = vertices[face[1]]
+                v2 = vertices[face[2]]
+                t = ray_intersect_triangle(origin, direction, v0, v1, v2)
+                if t is not None and t < best_t:
+                    best_t = t
+                    best_hit = origin + direction * t
+            if best_hit is not None:
+                points[y, x] = best_hit
+
+    return points, colors
+
+def save_colored_pointcloud(points, colors, output_path):
+    """Save as PLY file (colored point cloud)."""
+    # Flatten valid points
+    mask = ~np.isnan(points[:,:,0])
+    valid_points = points[mask].reshape(-1, 3)
+    valid_colors = (colors[mask].reshape(-1, 3) * 255).astype(np.uint8)
+
+    with open(output_path, 'w') as f:
+        f.write("ply\n")
+        f.write("format ascii 1.0\n")
+        f.write(f"element vertex {len(valid_points)}\n")
+        f.write("property float x\n")
+        f.write("property float y\n")
+        f.write("property float z\n")
+        f.write("property uchar red\n")
+        f.write("property uchar green\n")
+        f.write("property uchar blue\n")
+        f.write("end_header\n")
+        for pt, col in zip(valid_points, valid_colors):
+            f.write(f"{pt[0]} {pt[1]} {pt[2]} {col[0]} {col[1]} {col[2]}\n")
+    print(f"Saved point cloud to {output_path}")
+
+if __name__ == "__main__":
+    data_dir = "./exported_data"      # Path to exported files
+    image_path = "captured_image.png" # Input camera image
+    output_ply = "mapped_points.ply"
+
+    data = load_exported_data(data_dir)
+    points, colors = map_image_to_3d(image_path, data, use_precomputed=True)
+    save_colored_pointcloud(points, colors, output_ply)
+
+    # Optional: Visualize using Open3D
+    import open3d as o3d
+    pcd = o3d.geometry.PointCloud()
+    mask = ~np.isnan(points[:,:,0])
+    pcd.points = o3d.utility.Vector3dVector(points[mask].reshape(-1,3))
+    pcd.colors = o3d.utility.Vector3dVector(colors[mask].reshape(-1,3))
+    o3d.visualization.draw_geometries([pcd])

tools/3d_patch_generation_blender/test_blender.py

@@ -0,0 +1,278 @@
+import bpy
+import bmesh
+import numpy as np
+from mathutils import Vector
+import os
+
+# ----------------------------------------------------------------------
+# Scanline polygon fill for Blender Image pixels
+# ----------------------------------------------------------------------
+def fill_polygon(image, pixels, width, height, polygon, color):
+    poly = [(int(round(x)), int(round(y))) for x, y in polygon]
+    if len(poly) < 3:
+        return
+
+    min_x = max(0, min(p[0] for p in poly))
+    max_x = min(width-1, max(p[0] for p in poly))
+    min_y = max(0, min(p[1] for p in poly))
+    max_y = min(height-1, max(p[1] for p in poly))
+
+    for y in range(min_y, max_y + 1):
+        intersections = []
+        for i in range(len(poly)):
+            x1, y1 = poly[i]
+            x2, y2 = poly[(i+1) % len(poly)]
+            if y1 == y2:
+                continue
+            if (y1 <= y < y2) or (y2 <= y < y1):
+                t = (y - y1) / (y2 - y1)
+                x = x1 + t * (x2 - x1)
+                intersections.append(x)
+        intersections.sort()
+        for k in range(0, len(intersections), 2):
+            if k+1 >= len(intersections):
+                break
+            x_start = int(round(intersections[k]))
+            x_end = int(round(intersections[k+1]))
+            x_start = max(0, min(width-1, x_start))
+            x_end = max(0, min(width-1, x_end))
+            for x in range(x_start, x_end + 1):
+                idx = (y * width + x) * 4
+                pixels[idx:idx+4] = color
+
+# ----------------------------------------------------------------------
+# Get interpolated UV from face and hit location (fixed version)
+# ----------------------------------------------------------------------
+def get_uv_from_face(obj, face_index, location):
+    me = obj.data
+    local_loc = obj.matrix_world.inverted() @ location
+
+    poly = me.polygons[face_index]
+    vert_indices = poly.vertices
+    loop_indices = poly.loop_indices
+
+    if len(vert_indices) < 3:
+        return None
+
+    v0 = me.vertices[vert_indices[0]].co
+    v1 = me.vertices[vert_indices[1]].co
+    v2 = me.vertices[vert_indices[2]].co
+
+    uv_layer = me.uv_layers.active.data
+    uv0 = uv_layer[loop_indices[0]].uv
+    uv1 = uv_layer[loop_indices[1]].uv
+    uv2 = uv_layer[loop_indices[2]].uv
+
+    def tri_area(a, b, c):
+        return 0.5 * ((b - a).cross(c - a)).length
+
+    total_area = tri_area(v0, v1, v2)
+    if total_area == 0:
+        return None
+
+    w0 = tri_area(local_loc, v1, v2) / total_area
+    w1 = tri_area(local_loc, v2, v0) / total_area
+    w2 = 1.0 - w0 - w1
+
+    uv = w0 * uv0 + w1 * uv1 + w2 * uv2
+    return Vector((uv.x, uv.y))
+
+# ----------------------------------------------------------------------
+# Color generation functions (unchanged)
+# ----------------------------------------------------------------------
+def merge_intervals(intervals):
+    if not intervals:
+        return []
+    intervals.sort(key=lambda x: x[0])
+    merged = [intervals[0]]
+    for current in intervals[1:]:
+        last = merged[-1]
+        if current[0] <= last[1]:
+            merged[-1] = (last[0], max(last[1], current[1]))
+        else:
+            merged.append(current)
+    return merged
+
+def get_allowed_intervals(neighbor_values, delta):
+    forbidden = []
+    for c in neighbor_values:
+        low = max(0, c - delta)
+        high = min(1, c + delta)
+        forbidden.append((low, high))
+    merged = merge_intervals(forbidden)
+    allowed = []
+    if merged:
+        if merged[0][0] > 0:
+            allowed.append((0, merged[0][0]))
+        for i in range(1, len(merged)):
+            allowed.append((merged[i-1][1], merged[i][0]))
+        if merged[-1][1] < 1:
+            allowed.append((merged[-1][1], 1))
+    else:
+        allowed.append((0, 1))
+    return allowed
+
+def sample_from_allowed(allowed):
+    lengths = [end - start for start, end in allowed]
+    total_length = sum(lengths)
+    if total_length <= 0:
+        return None
+    point = np.random.rand() * total_length
+    cumulative = 0
+    for (start, end), length in zip(allowed, lengths):
+        cumulative += length
+        if point <= cumulative:
+            return np.random.uniform(start, end)
+    return allowed[-1][0]
+
+def sample_channel(neighbor_values, delta):
+    allowed = get_allowed_intervals(neighbor_values, delta)
+    if not allowed:
+        return None
+    return sample_from_allowed(allowed)
+
+def generate_valid_color(neighbors, r):
+    delta = np.sqrt(r)
+    color = np.zeros(3)
+    for channel in range(3):
+        neighbor_channel = [n[channel] for n in neighbors]
+        sampled = sample_channel(neighbor_channel, delta)
+        if sampled is None:
+            sampled = np.random.rand()
+        color[channel] = sampled
+    return color
+
+# ----------------------------------------------------------------------
+# Main function
+# ----------------------------------------------------------------------
+def generate_fingerprint_texture(
+    obj_name=None,
+    camera_name="Camera",
+    nx=40,
+    ny=30,
+    r=0.05,
+    texture_width=2048,
+    texture_height=2048,
+    output_path="//finger_texture.png"
+):
+    scene = bpy.context.scene
+    depsgraph = bpy.context.evaluated_depsgraph_get()
+    
+    if obj_name:
+        obj = bpy.data.objects[obj_name]
+    else:
+        obj = bpy.context.active_object
+    if obj is None or obj.type != 'MESH':
+        raise ValueError("No valid mesh object selected")
+    
+    cam = bpy.data.objects[camera_name]
+    if cam is None or cam.type != 'CAMERA':
+        raise ValueError("Camera not found")
+    
+    if not obj.data.uv_layers.active:
+        obj.data.uv_layers.new(name="UVMap")
+    
+    render = scene.render
+    res_x = render.resolution_x
+    res_y = render.resolution_y
+    
+    cam_data = cam.data
+    if cam_data.sensor_fit == 'HORIZONTAL':
+        sensor_width = cam_data.sensor_width
+        sensor_height = sensor_width / (res_x / res_y)
+    else:
+        sensor_height = cam_data.sensor_height
+        sensor_width = sensor_height * (res_x / res_y)
+    focal = cam_data.lens
+    frame_width = sensor_width / focal
+    frame_height = sensor_height / focal
+    
+    patch_uv_polygons = []
+    patch_colors = np.zeros((ny, nx, 3), dtype=np.float32)
+    
+    for j in range(ny):
+        v_start = j / ny
+        v_end = (j+1) / ny
+        for i in range(nx):
+            u_start = i / nx
+            u_end = (i+1) / nx
+            
+            corners_ndc = [
+                (u_start, v_start),
+                (u_end,   v_start),
+                (u_end,   v_end),
+                (u_start, v_end)
+            ]
+            
+            uv_poly = []
+            for (u_ndc, v_ndc) in corners_ndc:
+                sx = (u_ndc - 0.5) * frame_width
+                sy = (v_ndc - 0.5) * frame_height
+                direction_cam = Vector((sx, sy, -1.0))
+                direction_world = cam.matrix_world.to_quaternion() @ direction_cam
+                direction_world.normalize()
+                origin = cam.matrix_world.translation
+                
+                hit, location, normal, face_index, hit_obj, matrix = scene.ray_cast(depsgraph, origin, direction_world)
+                if hit and hit_obj == obj:
+                    uv = get_uv_from_face(obj, face_index, location)
+                    if uv is not None:
+                        uv_poly.append(uv)
+            
+            if len(uv_poly) < 3:
+                patch_uv_polygons.append(None)
+                continue
+            patch_uv_polygons.append(uv_poly)
+            
+            neighbors = []
+            if j > 0:
+                neighbors.append(patch_colors[j-1, i])
+            if i > 0:
+                neighbors.append(patch_colors[j, i-1])
+            
+            if len(neighbors) == 0:
+                color = np.random.rand(3)
+            else:
+                color = generate_valid_color(neighbors, r)
+            patch_colors[j, i] = color
+    
+    img_name = "FingerprintTexture"
+    if img_name in bpy.data.images:
+        bpy.data.images.remove(bpy.data.images[img_name])
+    image = bpy.data.images.new(img_name, width=texture_width, height=texture_height, alpha=True)
+    pixels = list(image.pixels)
+    
+    idx = 0
+    for j in range(ny):
+        for i in range(nx):
+            uv_poly = patch_uv_polygons[idx]
+            idx += 1
+            if uv_poly is None:
+                continue
+            pixel_poly = [(u * (texture_width - 1), v * (texture_height - 1)) for (u, v) in uv_poly]
+            color = patch_colors[j, i]
+            rgba = (color[0], color[1], color[2], 1.0)
+            fill_polygon(image, pixels, texture_width, texture_height, pixel_poly, rgba)
+    
+    image.pixels[:] = pixels
+    image.update()
+    
+    abs_path = bpy.path.abspath(output_path)
+    image.filepath_raw = abs_path
+    image.file_format = 'PNG'
+    image.save()
+    print(f"Texture saved to {abs_path}")
+    return image
+
+# ----------------------------------------------------------------------
+# Run
+# ----------------------------------------------------------------------
+if __name__ == "__main__":
+    generate_fingerprint_texture(
+        obj_name="FingerMesh",   # Change to your mesh name
+        nx=50, ny=40,
+        r=0.05,
+        texture_width=4096,
+        texture_height=4096,
+        output_path="//finger_texture.png"
+    )

tools/3d_patch_generation_blender/uv_to_printer.py

@@ -0,0 +1,180 @@
+"""
+Convert a UV texture map to a physically scaled image for UV printing.
+Usage: python uv_to_printer.py --model model.obj --texture texture.png --resolution 100 --output print_ready.png
+"""
+
+import argparse
+import numpy as np
+import cv2
+import trimesh
+import os
+
+def compute_uv_physical_scale(mesh, uv_layer_name=None):
+    """
+    Compute the physical dimensions (mm) of the UV space.
+    Assumes the UV map is a reasonably isometric unwrapping (e.g., cylindrical).
+    Returns: (width_mm, height_mm, scale_factor_u, scale_factor_v)
+        where scale_factor = physical length per UV unit.
+    """
+    # Get UV coordinates and corresponding 3D positions
+    if not hasattr(mesh.visual, 'uv') or mesh.visual.uv is None:
+        raise ValueError("Mesh has no UV coordinates.")
+
+    uvs = mesh.visual.uv
+    faces = mesh.faces
+    vertices = mesh.vertices
+
+    # For each face, we have 3 UVs and 3 vertex positions
+    # We'll compute the physical edge lengths for all edges in UV space
+    # and compare to 3D edge lengths to estimate scaling.
+
+    total_uv_length = 0.0
+    total_3d_length = 0.0
+    edge_count = 0
+
+    # Also track UV bounding box for overall dimensions
+    uv_min = uvs.min(axis=0)
+    uv_max = uvs.max(axis=0)
+    uv_range = uv_max - uv_min
+
+    # Estimate scale factor by sampling edges
+    for face in faces:
+        for i in range(3):
+            j = (i + 1) % 3
+            # UV edge
+            uv0 = uvs[face[i]]
+            uv1 = uvs[face[j]]
+            uv_len = np.linalg.norm(uv1 - uv0)
+            # 3D edge
+            v0 = vertices[face[i]]
+            v1 = vertices[face[j]]
+            len3d = np.linalg.norm(v1 - v0)
+
+            if uv_len > 1e-6:
+                total_uv_length += uv_len
+                total_3d_length += len3d
+                edge_count += 1
+
+    if edge_count == 0:
+        raise ValueError("No valid edges found for scale estimation.")
+
+    # Average scale factor (physical mm per UV unit)
+    scale_factor = total_3d_length / total_uv_length
+
+    # Physical dimensions of the UV layout
+    width_mm = uv_range[0] * scale_factor
+    height_mm = uv_range[1] * scale_factor
+
+    print(f"UV bounding box: {uv_range[0]:.3f} x {uv_range[1]:.3f} UV units")
+    print(f"Scale factor: {scale_factor:.3f} mm per UV unit")
+    print(f"Physical size: {width_mm:.2f} mm x {height_mm:.2f} mm")
+
+    return width_mm, height_mm, scale_factor
+
+def generate_print_image(texture_path, width_mm, height_mm, pixels_per_mm, output_path):
+    """
+    Load the texture image, scale it to match the physical dimensions at the
+    specified resolution, and save the result.
+    """
+    # Read texture (with alpha if present)
+    texture = cv2.imread(texture_path, cv2.IMREAD_UNCHANGED)
+    if texture is None:
+        raise FileNotFoundError(f"Cannot load texture: {texture_path}")
+
+    # Compute target pixel dimensions
+    target_width_px = int(round(width_mm * pixels_per_mm))
+    target_height_px = int(round(height_mm * pixels_per_mm))
+
+    print(f"Target image size: {target_width_px} x {target_height_px} pixels")
+
+    # Resize texture (using INTER_LINEAR for up/down scaling)
+    resized = cv2.resize(texture, (target_width_px, target_height_px),
+                         interpolation=cv2.INTER_LINEAR)
+
+    # Save output (preserve alpha if present)
+    cv2.imwrite(output_path, resized)
+    print(f"Print-ready image saved to: {output_path}")
+
+def process_exported_npz(npz_path, texture_path, pixels_per_mm, output_path):
+    """
+    Alternative: load data exported from Blender script (vertices, faces, UVs).
+    """
+    data = np.load(npz_path)
+    vertices = data['vertices']
+    faces = data['faces']
+    uv_data = data['uv_data']
+    loop_to_vertex = data['loop_to_vertex']
+
+    # Reconstruct per-vertex UVs (average)
+    # For simplicity, we'll compute physical scale using face edges as in trimesh method.
+    # We need to build a temporary mesh-like structure.
+
+    # Count number of faces
+    num_faces = len(faces)
+    # Build face UVs array (3 per face)
+    # The uv_data is per-loop, we need to group by face.
+    # Assuming faces are triangles.
+    face_uvs = []
+    face_vertices = []
+    uv_idx = 0
+    for face in faces:
+        f_uvs = []
+        f_verts = []
+        for _ in range(len(face)):
+            f_uvs.append(uv_data[uv_idx])
+            v_idx = loop_to_vertex[uv_idx]
+            f_verts.append(vertices[v_idx])
+            uv_idx += 1
+        face_uvs.append(np.array(f_uvs))
+        face_vertices.append(np.array(f_verts))
+
+    total_uv_len = 0.0
+    total_3d_len = 0.0
+    count = 0
+    uv_min = np.min(uv_data, axis=0)
+    uv_max = np.max(uv_data, axis=0)
+    uv_range = uv_max - uv_min
+
+    for f_uv, f_v in zip(face_uvs, face_vertices):
+        for i in range(3):
+            j = (i+1)%3
+            uv_len = np.linalg.norm(f_uv[i] - f_uv[j])
+            v_len = np.linalg.norm(f_v[i] - f_v[j])
+            if uv_len > 1e-6:
+                total_uv_len += uv_len
+                total_3d_len += v_len
+                count += 1
+
+    scale = total_3d_len / total_uv_len if count>0 else 1.0
+    width_mm = uv_range[0] * scale
+    height_mm = uv_range[1] * scale
+
+    print(f"Physical size from .npz: {width_mm:.2f} x {height_mm:.2f} mm")
+    generate_print_image(texture_path, width_mm, height_mm, pixels_per_mm, output_path)
+
+def main():
+    parser = argparse.ArgumentParser(description="Convert UV texture to physically scaled print image.")
+    parser.add_argument("--model", help="3D model file (OBJ, PLY, STL) with UVs")
+    parser.add_argument("--npz", help="Alternatively, use exported .npz from Blender script")
+    parser.add_argument("--texture", required=True, help="Input texture image (PNG, JPG)")
+    parser.add_argument("--resolution", type=float, default=100.0,
+                        help="Pixels per millimeter (default: 100)")
+    parser.add_argument("--output", default="print_ready.png", help="Output image path")
+    args = parser.parse_args()
+
+    if not args.model and not args.npz:
+        parser.error("Either --model or --npz must be provided.")
+
+    if args.npz:
+        process_exported_npz(args.npz, args.texture, args.resolution, args.output)
+    else:
+        # Load mesh with trimesh
+        mesh = trimesh.load(args.model, force='mesh')
+        if not isinstance(mesh, trimesh.Trimesh):
+            raise ValueError("Loaded object is not a single mesh.")
+
+        width_mm, height_mm, scale = compute_uv_physical_scale(mesh)
+        generate_print_image(args.texture, width_mm, height_mm, args.resolution, args.output)
+
+if __name__ == "__main__":
+    main()