Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
1.3k views
in Technique[技术] by (71.8m points)

math - 3d line-intersection code not working properly

I created this piece of code to get the intersection of two 3d line-segments.

Unfortunately the result of this code is inaccurate, the intersection-point is not always on both lines.

I am confused and unsure what I'm doing wrong.

Here is my code:

--dir = direction
--p1,p2 = represents the line
function GetIntersection(dirStart, dirEnd, p1, p2)
   local s1_x, s1_y, s2_x, s2_y =  dirEnd.x - dirStart.x, dirEnd.z - dirStart.z, p2.x - p1.x, p2.z - p1.z
   local div = (-s2_x * s1_y) + (s1_x * s2_y)

   if div == 0 then return nil end
   local s = (-s1_y * (dirStart.x - p1.x) + s1_x * (dirStart.z - p1.z)) / div
   local t = ( s2_x * (dirStart.z - p1.z) - s2_y * (dirStart.x - p1.x)) / div

   if (s >= 0 and s <= 1 and t >= 0 and t <= 1) and (Vector(dirStart.x + (t * s1_x), 0, dirStart.z + (t * s1_y)) or nil) then
      local v = Vector(dirStart.x + (t * s1_x),0,dirStart.z + (t * s1_y))
      return v
   end
end
See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

This is example of Delphi code to find a distance between two skew lines in 3D. For your purposes it is necessary to check that result if small enough value (intersection does exist), check that s and t parameters are in range 0..1, then calculate point using parameter s

Math of this approach is described in 'the shortest line...' section of Paul Bourke page

VecDiff if vector difference function, Dot id scalar product function

function LineLineDistance(const L0, L1: TLine3D; var s, t: Double): Double;
var
  u: TPoint3D;
  a, b, c, d, e, det, invdet:Double;
begin
  u := VecDiff(L1.Base, L0.Base);
  a := Dot(L0.Direction, L0.Direction);
  b := Dot(L0.Direction, L1.Direction);
  c := Dot(L1.Direction, L1.Direction);
  d := Dot(L0.Direction, u);
  e := Dot(L1.Direction, u);
  det := a * c - b * b;
  if det < eps then   
    Result := -1
  else begin
    invdet := 1 / det;
    s := invdet * (b * e - c * d);
    t := invdet * (a * e - b * d);

    Result := Distance(PointAtParam(L0, s), PointAtParam(L1, t));
  end;
end;

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...