§ Уравнение сферы

Сегодня рассмотрю как вычислить пересечение шара и линии, заданной параметрически.
Известно, что уравнение сферы задается как x^2 + y^2 + z^2 = r^2 . Здесь r - радиус сферы. Также, уравнение линии задается как M = A + t \vec{AB} . Здесь A - начальная позиция, а вектор - это вектор направления прямой.
Чтобы найти точку пересечения, надо x, y, z вставить в уравнение шара:
(A_x + tAB_x)^2 + (A_y + tAB_y)^2 + (A_z + tAB_z)^2 = r^2

§ Раскрытие скобок

Теперь необходимо раскрыть скобки:
(A_x + tAB_x)^2 = A_x^2 + 2tA_xAB_x + t^2AB_x^2
(A_y + tAB_y)^2 = A_y^2 + 2tA_yAB_y + t^2AB_y^2
(A_z + tAB_z)^2 = A_z^2 + 2tA_zAB_z + t^2AB_z^2
Приведем общие множители, и получится:
t^2 (AB_x^2 + AB_y^2 + AB_z^2) + 2t(A_xAB_x + A_yAB_y + A_zAB_z) + (A_x^2 + A_y^2 + A_z^2 - r^2) = 0
Теперь известны следующие коэффициенты:
A = AB_x^2 + AB_y^2 + AB_z^2
B = 2(A_xAB_x + A_yAB_y + A_zAB_z)
C = A_x^2 + A_y^2 + A_z^2 - r^2
Коэффициент С можно вычислить только один раз на фрейм.
Итого, получается квадратное уравнение:
at^2 + bt + c = 0

§ Решение квадратного уравнения

Чтобы решить квадратное уравнение, необходимо сначала найти дискриминант:
D = b^2 - 4ac
В зависимости от того, какой будет дискриминант, такое будет количество решений:
  • D < 0, решений нет, линия не имеет точки пересечения со сферой
  • D = 0, решение только одно
  • D > 0, решений два
Если есть либо одно решение, либо два, то находятся они так:
t_{1,2} = \frac{-b \pm \sqrt{D}}{2a}

§ Код

P - точка камеры, org - центр сферы, r - радиус.
1void render_sphere(vec3 P, vec3 org, float r) {
2
3    vec3  A = {P.x - org.x, P.y - org.y, P.z - org.z};
4    float c = A.x*A.x + A.y*A.y + A.z*A.z - r*r;
5
6    for (int y = -100; y < 100; y++)
7    for (int x = -160; x < 160; x++) {
8
9        vec3 AB = { (float)x/160, (float)y/160, 1.0 };
10
11        float a = AB.x*AB.x + AB.y*AB.y + AB.z*AB.z;
12        float b = 2*(A.x*AB.x + A.y*AB.y + A.z*AB.z);
13        float D = b*b - 4*a*c;
14
15        pset(160 + x, 100 - y, D > 0 ? 0xffffff : 0);
16    }
17}