Intersección de Rayo y Plano en C++

Abdul Mateen 12 octubre 2023
  1. Usar operaciones vectoriales en C++
  2. Intersección de Rayo y Plano en C++
Intersección de Rayo y Plano en C++

En este tutorial, primero, obtendremos una guía completa sobre la intersección del plano del rayo en C++.

Primero, discutiremos las operaciones vectoriales con su implementación. A continuación, discutiremos el concepto de intersección del plano del rayo y su implementación en C++.

Usar operaciones vectoriales en C++

Aquí, cubriremos las operaciones vectoriales relacionadas y su implementación en C++. Necesitamos las siguientes funciones:

  1. Operador menos: Para calcular la distancia entre dos puntos.
  2. Producto escalar: una función para calcular el producto escalar entre dos vectores, lo que da como resultado un número real.
  3. Operador de multiplicación: Para calcular el producto vectorial de dos vectores.
  4. Operador de multiplicación con parámetro flotante: para calcular la multiplicación escalar entre el vector y el valor escalar.
  5. Función de magnitud: Para calcular la magnitud del vector.
  6. Función de normalización: Para calcular la normal de un vector.
  7. Finalmente, el operador de flujo se sobrecarga para mostrar el vector.

Aquí está la clase 3D vectorial completa. La función main comprueba/demuestra funciones miembro de la clase 3D vectorial.

#include <math.h>

#include <iostream>

using namespace std;

class Vec3 {
  float x;
  float y;
  float z;

 public:
  Vec3() { x, y, z = 0; }
  Vec3(float x, float y, float z) {
    this->x = x;
    this->y = y;
    this->z = z;
  }
  Vec3 &operator+=(const Vec3 &b) {
    x = x + b.x;
    y = y + b.y;
    z = z + b.z;
    return *this;
  }
  Vec3 &operator-=(const Vec3 &b) {
    x = x - b.x;
    y = y - b.y;
    z = z - b.z;
    return *this;
  }
  Vec3 operator+(const Vec3 &b) {
    Vec3 newV = *this;
    newV += b;
    return newV;
  }
  Vec3 operator-(const Vec3 &b) {
    Vec3 newV = *this;
    newV -= b;
    return newV;
  }
  Vec3 operator*(const Vec3 &b) {  // cross operator
    Vec3 newV;
    newV.x = y * b.z - z * b.y;
    newV.y = x * b.z - z * b.x;
    newV.z = x * b.y - y * b.x;
    return newV;
  }
  Vec3 &operator*=(const float s) {  // Dot equal operator
    x = x * s;
    y = y * s;
    z = z * s;
    return *this;
  }
  Vec3 operator*(const float s) {  // Dot equal operator
    Vec3 newV = *this;
    return newV *= s;
  }
  float mag() { return sqrt(pow(x, 2) + pow(y, 2) + pow(z, 2)); }
  Vec3 &normalize() {
    double mag = this->mag();
    x /= mag;
    y /= mag;
    z /= mag;
    return *this;
  }
  // Dot operator
  float dot(const Vec3 &b) { return x * b.x + y * b.y + z * b.z; }
  friend ostream &operator<<(ostream &out, const Vec3 &v) {
    out << "(" << v.x << ", " << v.y << ", " << v.z << ")\n";
    return out;
  }
};

int main() {
  Vec3 v1(-5, 7, 2);
  Vec3 v2(4, 12, 1);
  Vec3 v3(1, 1, 1);
  cout << "Vec1 = " << v1 << endl;
  v1 += v2;
  cout << "Result of adding Vec2 in Vec1 = " << v1 << endl;
  v1 -= v3;
  cout << "Result of subtracting Vec3 from vec1 = " << v1 << endl;
  v1 *= 5;  // dot operation
  cout << "Resultant after the scaling of Vec1 with a value of 5 = " << v1
       << endl;
  cout << "Magnitude for Vec1 = " << v1.mag() << endl;
  v1.normalize();
  cout << "Vec1 after normalization = " << v1 << endl;
  cout << "Dot product of Vec1 and Vec2 = " << v1 * v2 << endl;
  return 0;
}

Aquí está el resultado de la función main, que muestra una implementación precisa de la clase 3D vectorial.

Vec1 = (-5, 7, 2)

Result of adding Vec2 in Vec1 = (-1, 19, 3)

Result of subtracting Vec3 from vec1 = (-2, 18, 2)

Resultant after the scaling of Vec1 with a value of 5 = (-10, 90, 10)

Magnitude for Vec1 = 91.1043
Vec1 after normalization = (-0.109764, 0.987878, 0.109764)

Dot product of Vec1 and Vec2 = (-0.329293, -0.548821, -5.26868)

Intersección de Rayo y Plano en C++

De aquí en adelante, asumimos que el lector se siente bastante cómodo con los conceptos de operaciones vectoriales. Suponemos además que el lector tiene una idea básica de los aviones.

Si no se siente cómodo con estos conceptos, siga este enlace para obtener más información sobre operaciones vectoriales y planos.

Veamos algunos hechos y detalles matemáticos formales para encontrar una intersección de rayos en el plano.

Sabemos que el producto escalar de dos vectores ortogonales o perpendiculares siempre es 0. Supongamos que a y b son dos vectores perpendiculares u ortogonales, entonces a.b=0.

Ahora, considere el punto p0 en el plano que representa la distancia del plano al origen y el vector n, normal al plano. Podemos calcular el vector p restando cualquier punto del plano del punto p0.

El vector resultante se encuentra en el plano y es perpendicular al plano normal.

Este hecho nos da la ecuación:

(p-p0) . N = 0                                            (i)

Considere l0 y l como el punto de partida del rayo y la dirección del rayo, respectivamente. Podemos llegar al plano (significa intersección en el punto p) usando la forma paramétrica:

l0 + l * t = p                                            (ii)

Si el rayo no es paralelo al plano, entonces t veces en la dirección del rayo; el rayo interesará al avión. A continuación, podemos poner el valor de p de la ecuación (ii) en la ecuación (i), y obtenemos:

(l0 + l * t – p0) . n = 0

Queremos calcular t, que puede ayudarnos a calcular la posición del punto de intersección usando la ecuación paramétrica. Resolvamos la ecuación:

l * t . n + (l0 – p0) . n = 0
l * t . n = - (l0 – p0) . n

t = - (l0 – p0) . n / l . n

Necesitamos determinar el producto escalar de la línea y el plano normales que vienen en el denominador porque si el valor es 0 o cercano a 0, el resultado será infinito, lo que significa que el rayo y el plano son paralelos. Por lo tanto, comprobaremos el denominador y devolveremos false, lo que significa que no hay solución (es decir, el punto de intersección).

La siguiente función puede verificar si el rayo se cruza con el plano o no:

bool intersection plane(Vec3 &n, Vec3 &p0, Vec3 &lo, Vec3 &l, float &t) {
  // assuming vectors are all normalized
  float denom = n.dot(l);
  if (denom < 0.00005 && denom > -0.00005) / denom is near to 0 return false;
  Vec3 p010 = p0  l0;
  t = p010.dot(n);
  if t>=0)    return true;
  return false;
}

Usando la condición if, comprobamos si el denominador está cerca de 0, lo que significa que habrá infinitos resultados. El producto escalar del plano normal con la dirección del rayo da 0 si el rayo es paralelo al plano; por lo tanto, devolvemos false.

De lo contrario, calcularemos t y devolveremos true, lo que significa que el rayo interseca al plano. Podemos encontrar un punto de intersección usando t.

A continuación, tenemos a continuación el código completo para encontrar el punto de intersección de un rayo y un plano.

#include <math.h>

#include <iostream>

using namespace std;

class Vec3 {
  float x;
  float y;
  float z;

 public:
  Vec3() { x, y, z = 0; }
  Vec3(float x, float y, float z) {
    this->x = x;
    this->y = y;
    this->z = z;
  }
  Vec3 &operator+=(const Vec3 &b) {
    x = x + b.x;
    y = y + b.y;
    z = z + b.z;
    return *this;
  }
  Vec3 &operator+=(const float s) {
    x = x + s;
    y = y + s;
    z = z + s;
    return *this;
  }
  Vec3 operator+(const float s) {
    Vec3 newV = *this;
    return newV += s;
  }
  Vec3 &operator-=(const Vec3 &b) {
    x = x - b.x;
    y = y - b.y;
    z = z - b.z;
    return *this;
  }
  Vec3 operator+(const Vec3 &b) {
    Vec3 newV = *this;
    newV += b;
    return newV;
  }
  Vec3 operator-(const Vec3 &b) {
    Vec3 newV = *this;
    newV -= b;
    return newV;
  }
  Vec3 operator*(const Vec3 &b) {  // cross operator
    Vec3 newV;
    newV.x = y * b.z - z * b.y;
    newV.y = x * b.z - z * b.x;
    newV.z = x * b.y - y * b.x;
    return newV;
  }
  Vec3 &operator*(const float s) {  // Dot equal operator
    x = x * s;
    y = y * s;
    z = z * s;
    return *this;
  }
  float mag() { return sqrt(pow(x, 2) + pow(y, 2) + pow(z, 2)); }
  Vec3 &normalize() {
    double mag = this->mag();
    x /= mag;
    y /= mag;
    z /= mag;
  }
  // Dot operator
  float dot(const Vec3 &b) { return x * b.x + y * b.y + z * b.z; }
  friend ostream &operator<<(ostream &out, const Vec3 &v) {
    out << "(" << v.x << ", " << v.y << ", " << v.z << ")\n";
    return out;
  }
};

bool intersectPlane(Vec3 &n, Vec3 &p0, Vec3 &l0, Vec3 &l, float &t) {
  // considering vectors are normalized
  float denom = n.dot(l);  // dot product n.l
  if (denom > 1e-6) {
    Vec3 p0l0 = p0 - l0;
    t = n.dot(p0l0) / denom;
    return (t >= 0);
  }
  return false;
}

int main() {
  Vec3 n1(3, -9, 1);  // Normal
  Vec3 p01(-4, 2, 2);
  Vec3 l01(1, 1, 1);
  Vec3 l1(1, 2, 1);
  float t;
  n1.normalize();
  p01.normalize();
  l01.normalize();
  l1.normalize();
  if (intersectPlane(n1, p01, l01, l1, t))
    cout << "T:" << t << '\n';
  else
    cout << "Ray is not intersecting the plane\n";
  cout << "------------------------\n";
  Vec3 n2(2, 2, -2);  // Normal
  Vec3 p02(2, 4, 1);
  Vec3 l02(1, 1, 1);
  Vec3 l2(1, 2, 1);
  n2.normalize();
  p02.normalize();
  l02.normalize();
  l2.normalize();
  if (intersectPlane(n2, p02, l02, l2, t))
    cout << "T:" << t << '\n';
  else
    cout << "Ray is not intersecting the plane\n";
  Vec3 intersectionPoint = l02 + l2 * t;
  cout << intersectionPoint;
  return 0;
}

La salida de este código es:

Ray is not intersecting the plane
------------------------
T:0.629199
(0.83422, 1.09109, 0.83422)

Como puede ver, en nuestro primer punto de datos, el rayo es paralelo al plano. Por lo tanto, el denominador es 0.

Sin embargo, en el segundo conjunto de puntos de datos, el denominador es mayor que 0, por lo que podemos encontrar el punto de intersección.

Artículo relacionado - C++ Math