Ray :: struct {
    origin      : Vector3;
    direction   : Vector3;
}

ray_aabb_intersection :: (origin: Vector3, direction: Vector3, aabb: AABB) -> float {
    t1 := (aabb.min.x - origin.x) / direction.x;
    t2 := (aabb.max.x - origin.x) / direction.x;
    t3 := (aabb.min.y - origin.y) / direction.y;
    t4 := (aabb.max.y - origin.y) / direction.y;
    t5 := (aabb.min.z - origin.z) / direction.z;
    t6 := (aabb.max.z - origin.z) / direction.z;

    tmin := max(max(min(t1, t2), min(t3, t4)), min(t5, t6));
    tmax := min(min(max(t1, t2), max(t3, t4)), max(t5, t6));

    // if tmax < 0, ray (line) is intersecting AABB, but whole AABB is behing us
    if tmax < 0 {
        return -1;
    }

    // if tmin > tmax, ray doesn't intersect AABB
    if tmin > tmax {
        return -1;
    }

    if tmin < 0 {
        return tmax;
    }
    return tmin;
}

ray_sphere_intersect :: (ray: Ray, center: Vector3, radius: float) -> bool {
    t : float;
    oc := ray.origin - center;
    a := dot(ray.direction, ray.direction);
    b := 2.0 * dot(oc, ray.direction);
    c := dot(oc, oc) - radius * radius;
    discriminant := b * b - 4 * a * c;

    if discriminant > 0 {
        // Two intersection points, choose the closest one (smallest positive t)
        t1 := (-b - sqrt(discriminant)) / (2.0 * a);
        t2 := (-b + sqrt(discriminant)) / (2.0 * a);

        if t1 > 0 || t2 > 0 {
            t = ifx (t1 < t2) then t1 else t2;
            return true;
        }
    } else if discriminant == 0 {
        // One intersection point
        t = -b / (2.0 * a);
        return t > 0;
    }

    // No intersection
    return false;
}

ray_plane_intersection :: (ray: Ray, p: Vector3, n: Vector3) -> bool, float {
    denom := dot(n, ray.direction);
    if abs(denom) > 0.0001 {
        offset := p - ray.origin;
        t := dot(offset, n) / denom;
        return t >= 0.0, t;
    }

    return false, 0.0;
}

ray_triangle_intersect :: (origin: Vector3, dir: Vector3, v0: Vector3, v1: Vector3, v2: Vector3) -> bool, float, Vector3 { 
    EPSILON := 0.0000001;
    edge1 := v1 - v0;
    edge2 := v2 - v0;
    h := cross(dir, edge2);
    a := dot(edge1, h);

    if a > -EPSILON && a < EPSILON return false, 0.0, .{0,0,0};    // This ray is parallel to this triangle.

    f := 1.0/a;
    s := origin - v0;
    u := f * dot(s, h);

    if u < 0.0 || u > 1.0 return false, 0.0, .{0,0,0};

    q := cross(s, edge1);
    v := f * dot(dir, q);

    if v < 0.0 || u + v > 1.0 return false, 0.0, .{0,0,0};

    // At this stage we can compute t to find out where the intersection point is on the line.
    t := f * dot(edge2, q);
    if t > EPSILON { // ray intersection
        //outIntersectionPoint = rayOrigin + rayVector * t;
        normal := normalize(cross(edge1, edge2));
        return true, t, normal;
    }
    else // This means that there is a line intersection but not a ray intersection.
        return false, 0.0, .{0,0,0};
    //e0 := v1 - v0;
    //e1 := v2 - v0;
    //pvec := cross(dir, e1);

    //det := dot(pvec, e0);

    //if det > 0 {
    //    inv_det := 1.0 / det;
    //    tvec := origin - v0;
    //    u := inv_det * dot(tvec, pvec);

    //    if u < 0.0 || u > 1.0 {
    //        return false, 0.0;
    //    }

    //    qvec := cross(tvec, e0);
    //    v := inv_det * dot(qvec, dir);

    //    if v < 0.0 || u + v > 1.0 {
    //        return false, 0.0;
    //    }

    //    return true, dot(e1, qvec) * inv_det;
    //}

    //return false, 0.0;
} 

distance_from_ray_to_point :: (ray: Ray, p: Vector3) -> float {
    pb := p - ray.origin;
    t0 := dot(ray.direction, pb) / dot(ray.direction, ray.direction);
    if  t0 <= 0.0 {
        return length(pb);
    } else if t0 > 0.0 && t0 < 1.0 {
        return length(p - (ray.origin + t0 * ray.direction));
    } else {
        return length(p - (ray.origin + ray.direction));
    }
}

closest_point_on_ray :: (from: Vector3, dir: Vector3, p: Vector3) -> Vector3 {
    lhs := p - from;
    dot_p := dot(lhs, dir);
    return from + dir * dot_p;
}

closest_distance_between_rays :: (l1: Ray, l2: Ray) -> distance: float, t1: float, t2: float {
    u := l1.direction;
    v := l2.direction;
    w := l1.origin - l2.origin;
    a := dot(u, u);
    b := dot(u, v);
    c := dot(v, v);
    d := dot(u, w);
    e := dot(v, w);
    D := a * c - b * b;

    sc: float;
    tc: float;

    if D < 0.00001 {
        sc = 0.0;

        if b > c {
            tc = d / b;
        } else {
            tc = e / c;
        }
    } else {
        sc = ((b * e) - (c * d)) / D;
        tc = ((a * e) - (b * d)) / D;
    }

    P1 := u * sc;
    P2 := v * tc;
    dP := w + P1 - P2;

    return length(dP), sc, tc;
}

get_positions :: (e : *Entity, n : Node_Render_Data, m : Mesh, i0 : u32, i1 : u32, i2 : u32) -> Vector3, Vector3, Vector3 {
	p0 : Vector3;
	p1 : Vector3;
	p2 : Vector3;
	
    p0 = transform_point(n.transform.world_matrix, m.positions[i0]);
    p1 = transform_point(n.transform.world_matrix, m.positions[i1]);
    p2 = transform_point(n.transform.world_matrix, m.positions[i2]);

	return p0, p1, p2;
}

ray_entity_intersect :: (ray: Ray, e: *Entity) -> bool, float, Vector3 {
    has_hit := false;
    closest : float = 10000000;
    closest_normal : Vector3;

    for n: e.renderable.model.nodes {
        render_node := e.renderable.nodes[it_index];
        for handle: n.meshes {
            m := parray_get(*engine.renderer.meshes, handle);
			index := 0;


			if m.indices.count > 0 {
				while index < m.indices.count - 1 {
					i0 := m.indices[index];
					i1 := m.indices[index + 1];
					i2 := m.indices[index + 2];

					p2, p1, p0 := get_positions(e, render_node, m, i0, i1, i2);

					result, t, normal := ray_triangle_intersect(ray.origin, ray.direction, p0, p1, p2);
					if result && t < closest {
                        has_hit = true;
                        closest = t;
                        closest_normal = normal;
					}

					index += 3;
				}
			} else {
				while index < m.positions.count - 1 {
					p2, p1, p0 := get_positions(e, render_node, m, xx index, xx (index + 1), xx (index + 2));

					result, t, normal := ray_triangle_intersect(ray.origin, ray.direction, p0, p1, p2);
					if result && t < closest {
                        has_hit = true;
                        closest = t;
                        closest_normal = normal;
					}

					index += 3;
				}
			}
        }
    }

    return has_hit, closest, closest_normal;
}