我需要一个基本函数来求点到线段的最短距离。你可以随意用任何你想要的语言来编写解决方案;我可以把它翻译成什么我正在使用(Javascript)。
编辑:线段由两个端点定义。线段AB由两点A (x1,y1)和B (x2,y2)定义。我要求的是这条线段到点C (x3,y3)的距离。我的几何技能生疏了,所以我看到的例子让我很困惑,我很遗憾地承认。
我需要一个基本函数来求点到线段的最短距离。你可以随意用任何你想要的语言来编写解决方案;我可以把它翻译成什么我正在使用(Javascript)。
编辑:线段由两个端点定义。线段AB由两点A (x1,y1)和B (x2,y2)定义。我要求的是这条线段到点C (x3,y3)的距离。我的几何技能生疏了,所以我看到的例子让我很困惑,我很遗憾地承认。
当前回答
Matlab代码,内置“自检”,如果他们调用函数没有参数:
function r = distPointToLineSegment( xy0, xy1, xyP )
% r = distPointToLineSegment( xy0, xy1, xyP )
if( nargin < 3 )
selfTest();
r=0;
else
vx = xy0(1)-xyP(1);
vy = xy0(2)-xyP(2);
ux = xy1(1)-xy0(1);
uy = xy1(2)-xy0(2);
lenSqr= (ux*ux+uy*uy);
detP= -vx*ux + -vy*uy;
if( detP < 0 )
r = norm(xy0-xyP,2);
elseif( detP > lenSqr )
r = norm(xy1-xyP,2);
else
r = abs(ux*vy-uy*vx)/sqrt(lenSqr);
end
end
function selfTest()
%#ok<*NASGU>
disp(['invalid args, distPointToLineSegment running (recursive) self-test...']);
ptA = [1;1]; ptB = [-1;-1];
ptC = [1/2;1/2]; % on the line
ptD = [-2;-1.5]; % too far from line segment
ptE = [1/2;0]; % should be same as perpendicular distance to line
ptF = [1.5;1.5]; % along the A-B but outside of the segment
distCtoAB = distPointToLineSegment(ptA,ptB,ptC)
distDtoAB = distPointToLineSegment(ptA,ptB,ptD)
distEtoAB = distPointToLineSegment(ptA,ptB,ptE)
distFtoAB = distPointToLineSegment(ptA,ptB,ptF)
figure(1); clf;
circle = @(x, y, r, c) rectangle('Position', [x-r, y-r, 2*r, 2*r], ...
'Curvature', [1 1], 'EdgeColor', c);
plot([ptA(1) ptB(1)],[ptA(2) ptB(2)],'r-x'); hold on;
plot(ptC(1),ptC(2),'b+'); circle(ptC(1),ptC(2), 0.5e-1, 'b');
plot(ptD(1),ptD(2),'g+'); circle(ptD(1),ptD(2), distDtoAB, 'g');
plot(ptE(1),ptE(2),'k+'); circle(ptE(1),ptE(2), distEtoAB, 'k');
plot(ptF(1),ptF(2),'m+'); circle(ptF(1),ptF(2), distFtoAB, 'm');
hold off;
axis([-3 3 -3 3]); axis equal;
end
end
其他回答
c#版本
public static FP DistanceToLineSegment(FPVector3 a, FPVector3 b, FPVector3 point)
{
var d = b - a;
var s = d.SqrMagnitude;
var ds = d / s;
var lambda = FPVector3.Dot(point - a, ds);
var p = FPMath.Clamp01(lambda) * d;
return (a + p - point).Magnitude;
}
对于感兴趣的人,这里是Joshua的Javascript代码到Objective-C的简单转换:
- (double)distanceToPoint:(CGPoint)p fromLineSegmentBetween:(CGPoint)l1 and:(CGPoint)l2
{
double A = p.x - l1.x;
double B = p.y - l1.y;
double C = l2.x - l1.x;
double D = l2.y - l1.y;
double dot = A * C + B * D;
double len_sq = C * C + D * D;
double param = dot / len_sq;
double xx, yy;
if (param < 0 || (l1.x == l2.x && l1.y == l2.y)) {
xx = l1.x;
yy = l1.y;
}
else if (param > 1) {
xx = l2.x;
yy = l2.y;
}
else {
xx = l1.x + param * C;
yy = l1.y + param * D;
}
double dx = p.x - xx;
double dy = p.y - yy;
return sqrtf(dx * dx + dy * dy);
}
我需要这个解决方案与MKMapPoint一起工作,所以我将分享它,以防其他人需要它。只是一些小的改变,这将返回米为单位的距离:
- (double)distanceToPoint:(MKMapPoint)p fromLineSegmentBetween:(MKMapPoint)l1 and:(MKMapPoint)l2
{
double A = p.x - l1.x;
double B = p.y - l1.y;
double C = l2.x - l1.x;
double D = l2.y - l1.y;
double dot = A * C + B * D;
double len_sq = C * C + D * D;
double param = dot / len_sq;
double xx, yy;
if (param < 0 || (l1.x == l2.x && l1.y == l2.y)) {
xx = l1.x;
yy = l1.y;
}
else if (param > 1) {
xx = l2.x;
yy = l2.y;
}
else {
xx = l1.x + param * C;
yy = l1.y + param * D;
}
return MKMetersBetweenMapPoints(p, MKMapPointMake(xx, yy));
}
只是遇到了这个,我想我应该添加一个Lua实现。它假设点以表{x=xVal, y=yVal}给出,直线或线段由包含两个点的表给出(见下面的例子):
function distance( P1, P2 )
return math.sqrt((P1.x-P2.x)^2 + (P1.y-P2.y)^2)
end
-- Returns false if the point lies beyond the reaches of the segment
function distPointToSegment( line, P )
if line[1].x == line[2].x and line[1].y == line[2].y then
print("Error: Not a line!")
return false
end
local d = distance( line[1], line[2] )
local t = ((P.x - line[1].x)*(line[2].x - line[1].x) + (P.y - line[1].y)*(line[2].y - line[1].y))/(d^2)
local projection = {}
projection.x = line[1].x + t*(line[2].x-line[1].x)
projection.y = line[1].y + t*(line[2].y-line[1].y)
if t >= 0 and t <= 1 then -- within line segment?
return distance( projection, {x=P.x, y=P.y} )
else
return false
end
end
-- Returns value even if point is further down the line (outside segment)
function distPointToLine( line, P )
if line[1].x == line[2].x and line[1].y == line[2].y then
print("Error: Not a line!")
return false
end
local d = distance( line[1], line[2] )
local t = ((P.x - line[1].x)*(line[2].x - line[1].x) + (P.y - line[1].y)*(line[2].y - line[1].y))/(d^2)
local projection = {}
projection.x = line[1].x + t*(line[2].x-line[1].x)
projection.y = line[1].y + t*(line[2].y-line[1].y)
return distance( projection, {x=P.x, y=P.y} )
end
使用示例:
local P1 = {x = 0, y = 0}
local P2 = {x = 10, y = 10}
local line = { P1, P2 }
local P3 = {x = 7, y = 15}
print(distPointToLine( line, P3 )) -- prints 5.6568542494924
print(distPointToSegment( line, P3 )) -- prints false
您可以尝试PHP geo-math-php的库
composer require rkondratuk/geo-math-php:^1
例子:
<?php
use PhpGeoMath\Model\GeoSegment;
use PhpGeoMath\Model\Polar3dPoint;
$polarPoint1 = new Polar3dPoint(
40.758742779050706, -73.97855507715238, Polar3dPoint::EARTH_RADIUS_IN_METERS
);
$polarPoint2 = new Polar3dPoint(
40.74843388072615, -73.98566565776102, Polar3dPoint::EARTH_RADIUS_IN_METERS
);
$polarPoint3 = new Polar3dPoint(
40.74919365249446, -73.98133456388013, Polar3dPoint::EARTH_RADIUS_IN_METERS
);
$arcSegment = new GeoSegment($polarPoint1, $polarPoint2);
$nearestPolarPoint = $arcSegment->calcNearestPoint($polarPoint3);
// Shortest distance from point-3 to segment(point-1, point-2)
$geoDistance = $nearestPolarPoint->calcGeoDistanceToPoint($polarPoint3);
Grumdrig的c++ /JavaScript实现对我来说非常有用,所以我提供了我正在使用的Python直接端口。完整的代码在这里。
class Point(object):
def __init__(self, x, y):
self.x = float(x)
self.y = float(y)
def square(x):
return x * x
def distance_squared(v, w):
return square(v.x - w.x) + square(v.y - w.y)
def distance_point_segment_squared(p, v, w):
# Segment length squared, |w-v|^2
d2 = distance_squared(v, w)
if d2 == 0:
# v == w, return distance to v
return distance_squared(p, v)
# Consider the line extending the segment, parameterized as v + t (w - v).
# We find projection of point p onto the line.
# It falls where t = [(p-v) . (w-v)] / |w-v|^2
t = ((p.x - v.x) * (w.x - v.x) + (p.y - v.y) * (w.y - v.y)) / d2;
if t < 0:
# Beyond v end of the segment
return distance_squared(p, v)
elif t > 1.0:
# Beyond w end of the segment
return distance_squared(p, w)
else:
# Projection falls on the segment.
proj = Point(v.x + t * (w.x - v.x), v.y + t * (w.y - v.y))
# print proj.x, proj.y
return distance_squared(p, proj)