📐

速度ベクトルを何倍すると直線に接触するか求める方法

2023/09/19に公開

目的

位置 c にいる物体が速度ベクトル cd で進んでいるとき、あと何倍すれば直線 ab と接触するかを知りたい。

考え方

次の 青 / 赤 をイメージする。

  • 青 → 「直線 ab の正規化」と「斜辺 ac」の外積で対辺の長さが求まる
  • 赤 → 「斜辺 cd」と「直線 ab の正規化」の外積で対辺の長さが求まる

計算方法

a = V[463.99999999999994, 60.0]
b = V[720.0, 340.0]
c = V[384.0, 252.0]
d = V[480.79999999999995, 246.0]

位置ベクトルから位置を外す。

cd = d - c  # => (96.79999999999995, -6.0)
ab = b - a  # => (256.00000000000006, 280.0)
ac = c - a  # => (-79.99999999999994, 192.0)
ad = d - a  # => (16.80000000000001, 186.0)
ca = a - c  # => (79.99999999999994, -192.0)

青は c からの垂線なので c_length とし、赤は d からの垂線なので d_length とする。そして比率から t が求まる。

c_length = ab.normalize.cross(ac)  # => 188.59799127280945
d_length = cd.cross(ab.normalize)  # => 75.48980419908962
t = c_length / d_length            # => 2.49832402234637

これで cd を約 2.5 倍すれば直線 ab に接触することがわかった。

高速化

t を求めるだけなら正規化を省略できる。上の計算方法では対辺の長さを求めるために ab を正規化したが、割合 t を求めるだけであれば分母・分子にある正規化は打ち消せる。つまり、次のように書いても t は同じになる。

c_length = ab.cross(ac)  # => 71552.0
d_length = cd.cross(ab)  # => 28639.999999999985
t = c_length / d_length  # => 2.49832402234637

一行で書けば

t = ab.cross(ac) / cd.cross(ab)  # => 2.49832402234637

となる。

公式で求める

まだ求まっていない交点を e とする。そして、二つのベクトルが平行だと外積は 0 になるので、ab × ae = 0と定義する。ae は ac と ce を足したものなので、ae = ac + ceと定義する。続いて ce は cd を 1 としたときの割合 t で求まるので、ce = cd * tと定義する。それを一つ前の式に代入すると、ae = ac + cd * tとなり、さらに前の式に代入すると、ab × (ac + cd * t) = 0となり、ここから t = の式にする。分配法則[1]により ab を両方に掛けて移項すると、(ab × ac) + (ab × cd * t) = 0ab × cd * t = -ab × act = -ab × ac / ab × cdとなる。ついでに外積は右辺と左辺を入れ替えると符号が反転するので、t = ab × ac / cd × abとしてマイナス符号を消す。できあがった t = の式を実際に計算すると

t = ab.cross(ac) / cd.cross(ab)  # => 2.49832402234637

となり 青 / 赤 をイメージしてできた数式と一致する。

内積を使って求める

青と赤は内積を使っても求めることができる。直線 ab を三角形の地面と考えたとき、青と赤の線は対辺なので、外積で求めることができた。ここで頭を90度回転させ直線 ab の法線をイメージする。その法線を地面と考えれば、青と赤の線は底辺なので、内積で求めることができる。

n = -ab.perp.normalize   # => (0.7380288120022732, -0.6747691995449356)
c_length = ca.dot(n)     # => 188.59799127280945
d_length = cd.dot(n)     # => 75.48980419908962
t = c_length / d_length  # => 2.49832402234637

法線は画面の右に伸びていると都合がよいので反転している。そして最後にこれも実際の長さが不要な場合は正規化を省略できるのと、法線の符号を考慮する必要もなくなるので

t = ca.dot(ab.perp) / cd.dot(ab.perp)  # => 2.49832402234637

と書ける。値は外積を使ったときと同じになっている。

交点を求める

もし交点 e が必要なら cd を t 倍する。

e = c.lerp(d, t)  # => (625.8377653631285, 237.01005586592177)

lerp は線形補完のことで

c + (d - c) * t  # => (625.8377653631285, 237.01005586592177)

と等しい。

関連

https://zenn.dev/megeton/articles/513247e5af6a46

参照

コード
class App < Base
  def initialize
    super

    self.width, self.height = V[800, 400]

    a = window_wh * V[0.580, 0.150]  # => (463.99999999999994, 60.0)
    b = window_wh * V[0.900, 0.850]  # => (720.0, 340.0)
    c = window_wh * V[0.480, 0.630]  # => (384.0, 252.0)
    d = window_wh * V[0.601, 0.615]  # => (480.79999999999995, 246.0)

    @points.concat([a, b, c, d])

    @mode = 0
    @perp = 0
  end

  def button_down(id)
    super

    if id == Gosu::KB_Z
      @mode = @mode.next.modulo(2)
    end
    if id == Gosu::KB_X
      @perp = @perp.next.modulo(2)
    end
  end

  def draw
    super

    a, b, c, d = @points
    cd = d - c
    ab = b - a
    ac = c - a
    ad = d - a
    aa = a - a
    ca = a - c
    t = ab.cross(ac) / cd.cross(ab)
    c_length = ab.normalize_or_zero.cross(ac)
    d_length = cd.cross(ab.normalize_or_zero)
    e = c.lerp(d, t)
    n = -ab.perp.normalize

    vputs "a: #{a.round}"
    vputs "b: #{b.round}"
    vputs "c: #{c.round}"
    vputs "d: #{d.round}"
    vputs ""

    vputs "公式"
    vputs "ab×ae = 0"
    vputs "ab×(ac + cd*t) = 0"
    vputs "(ab × ac) + (ab × cd*t) = 0"
    vputs "ab × cd*t = -ab × ac"
    vputs "t = -ab×ac / ab×cd"
    vputs "t = ab×ac / cd×ab = #{t.round(2)}"
    vputs ""

    vputs "青 = cの垂線 = abの正規化×ac = #{c_length.round}"
    vputs "赤 = dの垂線 = cd×abの正規化 = #{d_length.round}"
    vputs "青 / 赤 = #{(c_length / d_length).round(2)}"
    vputs ""

    vector_draw(a, b, "a", "b")
    vector_draw(c, d, "c", "d")

    if @perp == 1
      line_draw(c - ab.perp, c + ab.perp)
      line_draw(d - ab.perp, d + ab.perp)

      vputs "内積を使う場合"
      vputs "青 = cの垂線 = ca・abの法線の正規化 = #{ca.dot(n).round}"
      vputs "赤 = dの垂線 = cd・abの法線の正規化 = #{cd.dot(n).round}"
      vputs ""
    end

    if @mode == 0
    end

    if @mode == 1
      line_draw(a, c)
      line_draw(a - ab, b + ab)

      bold_line_draw(c, c - (ab.perp.normalize * ab.normalize_or_zero.cross(ac)), color: Gosu::Color::BLUE)
      bold_line_draw(d, d - (ab.perp.normalize * ab.normalize_or_zero.cross(cd)), color: Gosu::Color::RED)
      line_draw(c, c + cd.project_onto(ab))

      point_draw(c.lerp(d, t))
      arrow_head(c, e, "e")
    end
  end

  def font_size_default
    16
  end

  show
end
脚注
  1. A × (B + C) = A × B + A × C ↩︎

Discussion