C++ での光線と平面の交点

Abdul Mateen 2023年10月12日
  1. C++ でベクトル演算を使用する
  2. C++ での光線と平面の交点
C++ での光線と平面の交点

このチュートリアルでは、まず、C++ での光線と平面の交差に関する完全なガイダンスを取得します。

まず、ベクトル演算とその実装について説明します。 次に、光線平面交差の概念と C++ での実装について説明します。

C++ でベクトル演算を使用する

ここでは、関連するベクトル演算とその C++ での実装について説明します。 次の関数が必要です。

  1. マイナス演算子: 2 点間の距離を計算します。
  2. Dot product: 2つのベクトル間の内積を計算する関数で、結果は実数になります。
  3. 乗算演算子: 2つのベクトルの外積を計算します。
  4. float パラメーターを使用した乗算演算子: ベクトルとスカラー値の間のスカラー乗算を計算します。
  5. 大きさ関数: ベクトルの大きさを計算します。
  6. 正規化機能: ベクトルの法線を計算します。
  7. 最後に、ベクトルを表示するためにストリーム演算子がオーバーロードされます。

完全なベクター 3D クラスを次に示します。 main 関数は、ベクター 3D クラスのメンバー関数をチェック/デモンストレーションします。

#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;
}

以下は main 関数の出力で、ベクター 3D クラスの正確な実装を示しています。

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)

C++ での光線と平面の交点

ここから先は、読者がベクトル演算の概念に慣れていることを前提としています。 さらに、読者は平面について基本的な考えを持っていると仮定します。

これらの概念に慣れていない場合は、このリンク に従って、ベクトル操作と平面について詳しく学んでください。

平面内の光線の交差を見つけるために、いくつかの正式な数学的事実と詳細を見てみましょう。

2つの直交ベクトルまたは直交ベクトルの内積は常に 0 であることがわかっています。ab が 2つの直交ベクトルまたは直交ベクトルである場合、a.b=0 とします。

ここで、原点から平面までの距離を表す平面上の点p0と、平面に垂直なベクトルnを考えます。 点 p0 から平面上の任意の点を引くことにより、ベクトル p を計算できます。

結果のベクトルは平面内にあり、平面の法線に垂直です。

この事実から次の式が得られます。

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

l0l をそれぞれ光線の開始点と光線の方向とします。 パラメトリック形式を使用して、平面 (ポイント p で交差することを意味します) に到達できます。

l0 + l * t = p                                            (ii)

光線が平面に平行でない場合、光線の方向に t 回かかります。 光線は平面に関心を持ちます。 次に、式 (ii)p の値を式 (i) に代入すると、次のようになります。

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

t を計算します。これは、パラメトリック方程式を使用して交点の位置を計算するのに役立ちます。 方程式を解いてみましょう:

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

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

値が 0 または 0 に近い場合、結果は無限大になり、光線と平面が平行になるため、分母に来る線と平面の法線のドット積を決定する必要があります。 したがって、分母を確認して false を返します。これは、解がないことを意味します (つまり、交点)。

次の関数は、光線が平面と交差しているかどうかを確認できます。

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;
}

if 条件を使用して、分母が 0 に近いかどうかを確認します。これは、無限の結果があることを意味します。 光線が平面に平行な場合、平面法線と光線方向の内積は 0 になります。 したがって、false を返します。

それ以外の場合は、t を計算して true を返します。これは、光線が平面と交差することを意味します。 t を使用して交点を見つけることができます。

次に、光線と平面の交点を見つけるための完全なコードを以下に示します。

#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;
}

このコードの出力は次のとおりです。

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

ご覧のとおり、最初のデータ ポイントでは、光線は平面に平行です。 したがって、分母は 0 です。

ただし、2 番目のデータ ポイントのセットでは、分母が 0 より大きいため、交点を見つけることができます。

関連記事 - C++ Math