射线与椭球体求交的无依赖清晰代码实现方法咨询

射线与椭球求交完全可以在不依赖任何第三方图形库的前提下实现，核心逻辑采用仿射变换转单位球求交的方案，代码结构清晰且计算精度稳定。

核心计算流程 #

整个求交过程分三步完成：

1. 坐标变换：将世界坐标系下的射线平移、缩放到椭球的局部坐标系，将任意轴对齐椭球转换为球心在原点、半径为1的单位球
2. 单位球求交：求解射线与单位球的交点，过滤射线反方向的无效交点
3. 逆变换：将局部坐标系下的交点、法向量通过逆变换映射回世界坐标系

原有实现的问题 #

之前编写的球体判别式计算代码存在三个明显缺陷：

• 仅返回判别式值，未求解交点参数t，也未过滤t<0的无效交点（直线与球相交不代表射线与球相交）
• 未加入浮点数精度容错，射线起点位于曲面表面时极易出现自相交误判
• 缺失椭球到球体的仿射变换步骤，仅能支持正球体求交，无法适配椭球场景

完整实现代码 #

代码仅依赖C++标准库基础头文件，无任何图形库依赖：

#include <cmath>

// 三维向量/点结构体，自带基础运算
struct Vec3 {
    double x, y, z;
    Vec3() : x(0.0), y(0.0), z(0.0) {}
    Vec3(double x, double y, double z) : x(x), y(y), z(z) {}
    Vec3 operator-(const Vec3& other) const {
        return Vec3(x - other.x, y - other.y, z - other.z);
    }
    Vec3 operator+(const Vec3& other) const {
        return Vec3(x + other.x, y + other.y, z + other.z);
    }
    Vec3 operator*(double scalar) const {
        return Vec3(x * scalar, y * scalar, z * scalar);
    }
    Vec3 operator/(double scalar) const {
        return Vec3(x / scalar, y / scalar, z / scalar);
    }
};

// 向量点积工具函数
inline double dot(const Vec3& a, const Vec3& b) {
    return a.x * b.x + a.y * b.y + a.z * b.z;
}

// 椭球定义：中心坐标+三个正交轴的半轴长（轴对齐场景）
struct Ellipsoid {
    Vec3 center;
    double ra, rb, rc;
    Ellipsoid(Vec3 c, double a, double b, double c) : center(c), ra(a), rb(b), rc(c) {}
};

// 射线定义：起点+方向向量（方向无需提前归一化）
struct Ray {
    Vec3 origin;
    Vec3 dir;
    Ray(Vec3 o, Vec3 d) : origin(o), dir(d) {}
};

// 求交结果
struct HitRes {
    bool is_hit;
    double t; // 交点对应射线参数：hit_point = origin + t * dir
    Vec3 hit_point;
    Vec3 normal; // 世界空间归一化法向量
    HitRes() : is_hit(false), t(INFINITY) {}
};

HitRes rayIntersectEllipsoid(const Ray& ray, const Ellipsoid& e) {
    HitRes res;
    constexpr double EPS = 1e-6; // 浮点数精度阈值

    // 变换到椭球局部空间，椭球转换为单位球
    Vec3 local_orig = ray.origin - e.center;
    Vec3 scaled_orig(
        local_orig.x / e.ra,
        local_orig.y / e.rb,
        local_orig.z / e.rc
    );
    Vec3 scaled_dir(
        ray.dir.x / e.ra,
        ray.dir.y / e.rb,
        ray.dir.z / e.rc
    );

    // 构建单位球求交的一元二次方程 At² + Bt + C = 0
    double A = dot(scaled_dir, scaled_dir);
    double B = 2 * dot(scaled_orig, scaled_dir);
    double C = dot(scaled_orig, scaled_orig) - 1.0;

    double delta = B * B - 4 * A * C;
    if (delta < 0) return res; // 无实根，不相交

    double sqrt_delta = sqrt(delta);
    double t1 = (-B - sqrt_delta) / (2 * A);
    double t2 = (-B + sqrt_delta) / (2 * A);

    // 筛选最近的有效交点（t必须大于精度阈值，避免自相交）
    double t_hit = INFINITY;
    if (t1 > EPS) t_hit = t1;
    if (t2 > EPS && t2 < t_hit) t_hit = t2;
    if (t_hit == INFINITY) return res; // 所有交点都在射线后方

    // 局部空间交点逆变换回世界空间
    Vec3 local_hit = scaled_orig + scaled_dir * t_hit;
    res.is_hit = true;
    res.t = t_hit;
    res.hit_point = Vec3(
        local_hit.x * e.ra + e.center.x,
        local_hit.y * e.rb + e.center.y,
        local_hit.z * e.rc + e.center.z
    );

    // 计算世界空间法向量：法向量用逆转置变换，归一化后输出
    Vec3 world_normal(
        local_hit.x / e.ra,
        local_hit.y / e.rb,
        local_hit.z / e.rc
    );
    double normal_len = sqrt(dot(world_normal, world_normal));
    res.normal = world_normal / normal_len;

    return res;
}

调用示例 #

int main() {
    // 定义射线：起点(0,0,-5)，方向沿z轴正方向
    Ray ray(Vec3(0,0,-5), Vec3(0,0,1));
    // 定义椭球：中心在原点，x半轴1，y半轴2，z半轴1
    Ellipsoid e(Vec3(0,0,0), 1, 2, 1);

    HitRes res = rayIntersectEllipsoid(ray, e);
    if (res.is_hit) {
        //  res.hit_point 为交点坐标，res.normal为法向量，可直接用于渲染计算
    }
    return 0;
}

注意事项 #

• 上述实现针对轴对齐椭球，如果需要支持旋转椭球，在坐标变换步骤额外乘以旋转矩阵的逆矩阵即可
• 精度阈值1e-6可根据业务场景的数值范围调整，值过小容易出现自相交，值过大会丢失近距离的正确交点
• 传入的射线方向向量无需提前归一化，缩放逻辑会自动适配方向长度

内容的提问来源于stack exchange，提问作者DumbSimon