Source code for simvx.core.math.raycast

"""Raycasting helpers: screen-to-ray and ray-sphere intersection."""

import math

import numpy as np

from .types import Vec3




def screen_to_ray(screen_pos, screen_size, view, proj):
    """Convert a screen pixel coordinate to a world-space ray (origin, direction).

    Returns (origin, direction) as Vec3.
    """
    sp_x = float(screen_pos[0]) if isinstance(screen_pos, tuple | list) else float(screen_pos.x)
    sp_y = float(screen_pos[1]) if isinstance(screen_pos, tuple | list) else float(screen_pos.y)
    ss_x = float(screen_size[0]) if isinstance(screen_size, tuple | list) else float(screen_size.x)
    ss_y = float(screen_size[1]) if isinstance(screen_size, tuple | list) else float(screen_size.y)
    ndc_x = (2.0 * sp_x / ss_x) - 1.0
    ndc_y = 1.0 - (2.0 * sp_y / ss_y)

    if not isinstance(view, np.ndarray):
        view = np.array(view, dtype=np.float32).reshape(4, 4)
    if not isinstance(proj, np.ndarray):
        proj = np.array(proj, dtype=np.float32).reshape(4, 4)

    # Detect Vulkan Y-flip
    if proj[1, 1] < 0:
        ndc_y = -ndc_y

    inv_vp = np.linalg.inv(proj @ view)
    near = inv_vp @ np.array([ndc_x, ndc_y, -1.0, 1.0], dtype=np.float32)
    far = inv_vp @ np.array([ndc_x, ndc_y, 1.0, 1.0], dtype=np.float32)
    near_pos = near[:3] / near[3]
    far_pos = far[:3] / far[3]
    direction = far_pos - near_pos
    direction = direction / np.linalg.norm(direction)
    return Vec3(near_pos), Vec3(direction)





def ray_intersect_sphere(origin, direction, center, radius: float):
    """Ray-sphere intersection. Returns distance t or None if no hit."""
    if not isinstance(origin, Vec3):
        origin = Vec3(origin)
    if not isinstance(direction, Vec3):
        direction = Vec3(direction)
    if not isinstance(center, Vec3):
        center = Vec3(center)
    oc = origin - center
    b = oc.dot(direction)
    c = oc.dot(oc) - radius * radius
    discriminant = b * b - c
    if discriminant < 0:
        return None
    sqrt_d = math.sqrt(discriminant)
    t = -b - sqrt_d
    if t < 0:
        t = -b + sqrt_d
    return t if t >= 0 else None