我有一条从a到B的直线和一个半径为R的圆。
用什么算法来检查直线是否与圆相交?它在圆边的哪个坐标上?
我有一条从a到B的直线和一个半径为R的圆。
用什么算法来检查直线是否与圆相交?它在圆边的哪个坐标上?
当前回答
我会用这个算法来计算点(圆心)和线(线AB)之间的距离。这可以用来确定直线与圆的交点。
假设有点A B c, Ax和Ay是A点的x和y分量。B和c也是一样,标量R是圆半径。
该算法要求A B C是不同的点,且R不为0。
这是算法
// compute the euclidean distance between A and B
LAB = sqrt( (Bx-Ax)²+(By-Ay)² )
// compute the direction vector D from A to B
Dx = (Bx-Ax)/LAB
Dy = (By-Ay)/LAB
// the equation of the line AB is x = Dx*t + Ax, y = Dy*t + Ay with 0 <= t <= LAB.
// compute the distance between the points A and E, where
// E is the point of AB closest the circle center (Cx, Cy)
t = Dx*(Cx-Ax) + Dy*(Cy-Ay)
// compute the coordinates of the point E
Ex = t*Dx+Ax
Ey = t*Dy+Ay
// compute the euclidean distance between E and C
LEC = sqrt((Ex-Cx)²+(Ey-Cy)²)
// test if the line intersects the circle
if( LEC < R )
{
// compute distance from t to circle intersection point
dt = sqrt( R² - LEC²)
// compute first intersection point
Fx = (t-dt)*Dx + Ax
Fy = (t-dt)*Dy + Ay
// compute second intersection point
Gx = (t+dt)*Dx + Ax
Gy = (t+dt)*Dy + Ay
}
// else test if the line is tangent to circle
else if( LEC == R )
// tangent point to circle is E
else
// line doesn't touch circle
其他回答
我根据chmike给出的答案为iOS创建了这个函数
+ (NSArray *)intersectionPointsOfCircleWithCenter:(CGPoint)center withRadius:(float)radius toLinePoint1:(CGPoint)p1 andLinePoint2:(CGPoint)p2
{
NSMutableArray *intersectionPoints = [NSMutableArray array];
float Ax = p1.x;
float Ay = p1.y;
float Bx = p2.x;
float By = p2.y;
float Cx = center.x;
float Cy = center.y;
float R = radius;
// compute the euclidean distance between A and B
float LAB = sqrt( pow(Bx-Ax, 2)+pow(By-Ay, 2) );
// compute the direction vector D from A to B
float Dx = (Bx-Ax)/LAB;
float Dy = (By-Ay)/LAB;
// Now the line equation is x = Dx*t + Ax, y = Dy*t + Ay with 0 <= t <= 1.
// compute the value t of the closest point to the circle center (Cx, Cy)
float t = Dx*(Cx-Ax) + Dy*(Cy-Ay);
// This is the projection of C on the line from A to B.
// compute the coordinates of the point E on line and closest to C
float Ex = t*Dx+Ax;
float Ey = t*Dy+Ay;
// compute the euclidean distance from E to C
float LEC = sqrt( pow(Ex-Cx, 2)+ pow(Ey-Cy, 2) );
// test if the line intersects the circle
if( LEC < R )
{
// compute distance from t to circle intersection point
float dt = sqrt( pow(R, 2) - pow(LEC,2) );
// compute first intersection point
float Fx = (t-dt)*Dx + Ax;
float Fy = (t-dt)*Dy + Ay;
// compute second intersection point
float Gx = (t+dt)*Dx + Ax;
float Gy = (t+dt)*Dy + Ay;
[intersectionPoints addObject:[NSValue valueWithCGPoint:CGPointMake(Fx, Fy)]];
[intersectionPoints addObject:[NSValue valueWithCGPoint:CGPointMake(Gx, Gy)]];
}
// else test if the line is tangent to circle
else if( LEC == R ) {
// tangent point to circle is E
[intersectionPoints addObject:[NSValue valueWithCGPoint:CGPointMake(Ex, Ey)]];
}
else {
// line doesn't touch circle
}
return intersectionPoints;
}
采取
E是射线的起点, L是射线的端点, C是你测试的圆心 R是球面的半径
计算: d = L - E(射线方向矢量,从头到尾) f = E - C(从中心球到射线起点的向量)
然后通过…找到交点。 堵塞: P = E + t * d 这是一个参数方程 Px = Ex + tdx Py = Ey + tdy 成 (x - h)2 + (y - k)2 = r2 (h,k) =圆心。
注意:我们在这里将问题简化为2D,我们得到的解决方案也适用于3D
得到:
Expand x2 - 2xh + h2 + y2 - 2yk + k2 - r2 = 0 Plug x = ex + tdx y = ey + tdy ( ex + tdx )2 - 2( ex + tdx )h + h2 + ( ey + tdy )2 - 2( ey + tdy )k + k2 - r2 = 0 Explode ex2 + 2extdx + t2dx2 - 2exh - 2tdxh + h2 + ey2 + 2eytdy + t2dy2 - 2eyk - 2tdyk + k2 - r2 = 0 Group t2( dx2 + dy2 ) + 2t( exdx + eydy - dxh - dyk ) + ex2 + ey2 - 2exh - 2eyk + h2 + k2 - r2 = 0 Finally, t2( d · d ) + 2t( e · d - d · c ) + e · e - 2( e · c ) + c · c - r2 = 0 Where d is the vector d and · is the dot product. And then, t2( d · d ) + 2t( d · ( e - c ) ) + ( e - c ) · ( e - c ) - r2 = 0 Letting f = e - c t2( d · d ) + 2t( d · f ) + f · f - r2 = 0
所以我们得到: T2 *(d·d) + 2t*(f·d) + (f·f - r2) = 0
求解二次方程:
float a = d.Dot( d ) ;
float b = 2*f.Dot( d ) ;
float c = f.Dot( f ) - r*r ;
float discriminant = b*b-4*a*c;
if( discriminant < 0 )
{
// no intersection
}
else
{
// ray didn't totally miss sphere,
// so there is a solution to
// the equation.
discriminant = sqrt( discriminant );
// either solution may be on or off the ray so need to test both
// t1 is always the smaller value, because BOTH discriminant and
// a are nonnegative.
float t1 = (-b - discriminant)/(2*a);
float t2 = (-b + discriminant)/(2*a);
// 3x HIT cases:
// -o-> --|--> | | --|->
// Impale(t1 hit,t2 hit), Poke(t1 hit,t2>1), ExitWound(t1<0, t2 hit),
// 3x MISS cases:
// -> o o -> | -> |
// FallShort (t1>1,t2>1), Past (t1<0,t2<0), CompletelyInside(t1<0, t2>1)
if( t1 >= 0 && t1 <= 1 )
{
// t1 is the intersection, and it's closer than t2
// (since t1 uses -b - discriminant)
// Impale, Poke
return true ;
}
// here t1 didn't intersect so we are either started
// inside the sphere or completely past it
if( t2 >= 0 && t2 <= 1 )
{
// ExitWound
return true ;
}
// no intn: FallShort, Past, CompletelyInside
return false ;
}
另一个在c#(部分圆类)。 经过测试,工作就像一个魅力。
public class Circle : IEquatable<Circle>
{
// ******************************************************************
// The center of a circle
private Point _center;
// The radius of a circle
private double _radius;
// ******************************************************************
/// <summary>
/// Find all intersections (0, 1, 2) of the circle with a line defined by its 2 points.
/// Using: http://math.stackexchange.com/questions/228841/how-do-i-calculate-the-intersections-of-a-straight-line-and-a-circle
/// Note: p is the Center.X and q is Center.Y
/// </summary>
/// <param name="linePoint1"></param>
/// <param name="linePoint2"></param>
/// <returns></returns>
public List<Point> GetIntersections(Point linePoint1, Point linePoint2)
{
List<Point> intersections = new List<Point>();
double dx = linePoint2.X - linePoint1.X;
if (dx.AboutEquals(0)) // Straight vertical line
{
if (linePoint1.X.AboutEquals(Center.X - Radius) || linePoint1.X.AboutEquals(Center.X + Radius))
{
Point pt = new Point(linePoint1.X, Center.Y);
intersections.Add(pt);
}
else if (linePoint1.X > Center.X - Radius && linePoint1.X < Center.X + Radius)
{
double x = linePoint1.X - Center.X;
Point pt = new Point(linePoint1.X, Center.Y + Math.Sqrt(Radius * Radius - (x * x)));
intersections.Add(pt);
pt = new Point(linePoint1.X, Center.Y - Math.Sqrt(Radius * Radius - (x * x)));
intersections.Add(pt);
}
return intersections;
}
// Line function (y = mx + b)
double dy = linePoint2.Y - linePoint1.Y;
double m = dy / dx;
double b = linePoint1.Y - m * linePoint1.X;
double A = m * m + 1;
double B = 2 * (m * b - m * _center.Y - Center.X);
double C = Center.X * Center.X + Center.Y * Center.Y - Radius * Radius - 2 * b * Center.Y + b * b;
double discriminant = B * B - 4 * A * C;
if (discriminant < 0)
{
return intersections; // there is no intersections
}
if (discriminant.AboutEquals(0)) // Tangeante (touch on 1 point only)
{
double x = -B / (2 * A);
double y = m * x + b;
intersections.Add(new Point(x, y));
}
else // Secant (touch on 2 points)
{
double x = (-B + Math.Sqrt(discriminant)) / (2 * A);
double y = m * x + b;
intersections.Add(new Point(x, y));
x = (-B - Math.Sqrt(discriminant)) / (2 * A);
y = m * x + b;
intersections.Add(new Point(x, y));
}
return intersections;
}
// ******************************************************************
// Get the center
[XmlElement("Center")]
public Point Center
{
get { return _center; }
set
{
_center = value;
}
}
// ******************************************************************
// Get the radius
[XmlElement]
public double Radius
{
get { return _radius; }
set { _radius = value; }
}
//// ******************************************************************
//[XmlArrayItemAttribute("DoublePoint")]
//public List<Point> Coordinates
//{
// get { return _coordinates; }
//}
// ******************************************************************
// Construct a circle without any specification
public Circle()
{
_center.X = 0;
_center.Y = 0;
_radius = 0;
}
// ******************************************************************
// Construct a circle without any specification
public Circle(double radius)
{
_center.X = 0;
_center.Y = 0;
_radius = radius;
}
// ******************************************************************
// Construct a circle with the specified circle
public Circle(Circle circle)
{
_center = circle._center;
_radius = circle._radius;
}
// ******************************************************************
// Construct a circle with the specified center and radius
public Circle(Point center, double radius)
{
_center = center;
_radius = radius;
}
// ******************************************************************
// Construct a circle based on one point
public Circle(Point center)
{
_center = center;
_radius = 0;
}
// ******************************************************************
// Construct a circle based on two points
public Circle(Point p1, Point p2)
{
Circle2Points(p1, p2);
}
要求:
using System;
namespace Mathematic
{
public static class DoubleExtension
{
// ******************************************************************
// Base on Hans Passant Answer on:
// http://stackoverflow.com/questions/2411392/double-epsilon-for-equality-greater-than-less-than-less-than-or-equal-to-gre
/// <summary>
/// Compare two double taking in account the double precision potential error.
/// Take care: truncation errors accumulate on calculation. More you do, more you should increase the epsilon.
public static bool AboutEquals(this double value1, double value2)
{
if (double.IsPositiveInfinity(value1))
return double.IsPositiveInfinity(value2);
if (double.IsNegativeInfinity(value1))
return double.IsNegativeInfinity(value2);
if (double.IsNaN(value1))
return double.IsNaN(value2);
double epsilon = Math.Max(Math.Abs(value1), Math.Abs(value2)) * 1E-15;
return Math.Abs(value1 - value2) <= epsilon;
}
// ******************************************************************
// Base on Hans Passant Answer on:
// http://stackoverflow.com/questions/2411392/double-epsilon-for-equality-greater-than-less-than-less-than-or-equal-to-gre
/// <summary>
/// Compare two double taking in account the double precision potential error.
/// Take care: truncation errors accumulate on calculation. More you do, more you should increase the epsilon.
/// You get really better performance when you can determine the contextual epsilon first.
/// </summary>
/// <param name="value1"></param>
/// <param name="value2"></param>
/// <param name="precalculatedContextualEpsilon"></param>
/// <returns></returns>
public static bool AboutEquals(this double value1, double value2, double precalculatedContextualEpsilon)
{
if (double.IsPositiveInfinity(value1))
return double.IsPositiveInfinity(value2);
if (double.IsNegativeInfinity(value1))
return double.IsNegativeInfinity(value2);
if (double.IsNaN(value1))
return double.IsNaN(value2);
return Math.Abs(value1 - value2) <= precalculatedContextualEpsilon;
}
// ******************************************************************
public static double GetContextualEpsilon(this double biggestPossibleContextualValue)
{
return biggestPossibleContextualValue * 1E-15;
}
// ******************************************************************
/// <summary>
/// Mathlab equivalent
/// </summary>
/// <param name="dividend"></param>
/// <param name="divisor"></param>
/// <returns></returns>
public static double Mod(this double dividend, double divisor)
{
return dividend - System.Math.Floor(dividend / divisor) * divisor;
}
// ******************************************************************
}
}
You can find a point on a infinite line that is nearest to circle center by projecting vector AC onto vector AB. Calculate the distance between that point and circle center. If it is greater that R, there is no intersection. If the distance is equal to R, line is a tangent of the circle and the point nearest to circle center is actually the intersection point. If distance less that R, then there are 2 intersection points. They lie at the same distance from the point nearest to circle center. That distance can easily be calculated using Pythagorean theorem. Here's algorithm in pseudocode:
{
dX = bX - aX;
dY = bY - aY;
if ((dX == 0) && (dY == 0))
{
// A and B are the same points, no way to calculate intersection
return;
}
dl = (dX * dX + dY * dY);
t = ((cX - aX) * dX + (cY - aY) * dY) / dl;
// point on a line nearest to circle center
nearestX = aX + t * dX;
nearestY = aY + t * dY;
dist = point_dist(nearestX, nearestY, cX, cY);
if (dist == R)
{
// line segment touches circle; one intersection point
iX = nearestX;
iY = nearestY;
if (t < 0 || t > 1)
{
// intersection point is not actually within line segment
}
}
else if (dist < R)
{
// two possible intersection points
dt = sqrt(R * R - dist * dist) / sqrt(dl);
// intersection point nearest to A
t1 = t - dt;
i1X = aX + t1 * dX;
i1Y = aY + t1 * dY;
if (t1 < 0 || t1 > 1)
{
// intersection point is not actually within line segment
}
// intersection point farthest from A
t2 = t + dt;
i2X = aX + t2 * dX;
i2Y = aY + t2 * dY;
if (t2 < 0 || t2 > 1)
{
// intersection point is not actually within line segment
}
}
else
{
// no intersection
}
}
编辑:增加了代码来检查所找到的交点是否实际上在线段内。
似乎没人考虑投影,我是不是完全跑题了?
将向量AC投影到AB上,投影的向量AD就得到了新的点D。 如果D和C之间的距离小于(或等于)R,我们有一个交点。
是这样的:
社区编辑:
对于稍后无意中看到这篇文章并想知道如何实现这样一个算法的人来说,这里是一个使用常见向量操作函数用JavaScript编写的通用实现。
/**
* Returns the distance from line segment AB to point C
*/
function distanceSegmentToPoint(A, B, C) {
// Compute vectors AC and AB
const AC = sub(C, A);
const AB = sub(B, A);
// Get point D by taking the projection of AC onto AB then adding the offset of A
const D = add(proj(AC, AB), A);
const AD = sub(D, A);
// D might not be on AB so calculate k of D down AB (aka solve AD = k * AB)
// We can use either component, but choose larger value to reduce the chance of dividing by zero
const k = Math.abs(AB.x) > Math.abs(AB.y) ? AD.x / AB.x : AD.y / AB.y;
// Check if D is off either end of the line segment
if (k <= 0.0) {
return Math.sqrt(hypot2(C, A));
} else if (k >= 1.0) {
return Math.sqrt(hypot2(C, B));
}
return Math.sqrt(hypot2(C, D));
}
对于这个实现,我使用了两个常见的矢量操作函数,无论您在什么环境中工作,都可能已经提供了这些函数。但是,如果您还没有这些可用的功能,下面介绍如何实现它们。
// Define some common functions for working with vectors
const add = (a, b) => ({x: a.x + b.x, y: a.y + b.y});
const sub = (a, b) => ({x: a.x - b.x, y: a.y - b.y});
const dot = (a, b) => a.x * b.x + a.y * b.y;
const hypot2 = (a, b) => dot(sub(a, b), sub(a, b));
// Function for projecting some vector a onto b
function proj(a, b) {
const k = dot(a, b) / dot(b, b);
return {x: k * b.x, y: k * b.y};
}