PostGIS  3.0.0dev-r@@SVN_REVISION@@

◆ lw_dist2d_arc_arc()

int lw_dist2d_arc_arc ( const POINT2D A1,
const POINT2D A2,
const POINT2D A3,
const POINT2D B1,
const POINT2D B2,
const POINT2D B3,
DISTPTS dl 
)

Definition at line 1497 of file measures.c.

References DIST_MIN, DISTPTS::distance, distance2d_pt_pt(), FP_EQUALS, lw_arc_center(), lw_arc_is_pt(), lw_dist2d_arc_arc_concentric(), lw_dist2d_pt_arc(), lw_dist2d_pt_pt(), lw_dist2d_seg_arc(), lw_dist2d_seg_seg(), LW_FALSE, lw_pt_in_arc(), LW_TRUE, lwerror(), DISTPTS::mode, DISTPTS::p1, DISTPTS::p2, POINT2D::x, and POINT2D::y.

Referenced by lw_dist2d_ptarrayarc_ptarrayarc(), rect_leaf_node_distance(), rect_leaf_node_intersects(), and test_lw_dist2d_arc_arc().

1500 {
1501  POINT2D CA, CB; /* Center points of arcs A and B */
1502  double radius_A, radius_B, d; /* Radii of arcs A and B */
1503  POINT2D P; /* Temporary point P */
1504  POINT2D D; /* Mid-point between the centers CA and CB */
1505  int pt_in_arc_A, pt_in_arc_B; /* Test whether potential intersection point is within the arc */
1506 
1507  if ( dl->mode != DIST_MIN )
1508  lwerror("lw_dist2d_arc_arc only supports mindistance");
1509 
1510  /* TODO: Handle case where arc is closed circle (A1 = A3) */
1511 
1512  /* What if one or both of our "arcs" is actually a point? */
1513  if ( lw_arc_is_pt(B1, B2, B3) && lw_arc_is_pt(A1, A2, A3) )
1514  return lw_dist2d_pt_pt(B1, A1, dl);
1515  else if ( lw_arc_is_pt(B1, B2, B3) )
1516  return lw_dist2d_pt_arc(B1, A1, A2, A3, dl);
1517  else if ( lw_arc_is_pt(A1, A2, A3) )
1518  return lw_dist2d_pt_arc(A1, B1, B2, B3, dl);
1519 
1520  /* Calculate centers and radii of circles. */
1521  radius_A = lw_arc_center(A1, A2, A3, &CA);
1522  radius_B = lw_arc_center(B1, B2, B3, &CB);
1523 
1524  /* Two co-linear arcs?!? That's two segments. */
1525  if ( radius_A < 0 && radius_B < 0 )
1526  return lw_dist2d_seg_seg(A1, A3, B1, B3, dl);
1527 
1528  /* A is co-linear, delegate to lw_dist_seg_arc here. */
1529  if ( radius_A < 0 )
1530  return lw_dist2d_seg_arc(A1, A3, B1, B2, B3, dl);
1531 
1532  /* B is co-linear, delegate to lw_dist_seg_arc here. */
1533  if ( radius_B < 0 )
1534  return lw_dist2d_seg_arc(B1, B3, A1, A2, A3, dl);
1535 
1536  /* Center-center distance */
1537  d = distance2d_pt_pt(&CA, &CB);
1538 
1539  /* Concentric arcs */
1540  if ( FP_EQUALS(d, 0.0) )
1541  return lw_dist2d_arc_arc_concentric(A1, A2, A3, radius_A,
1542  B1, B2, B3, radius_B,
1543  &CA, dl);
1544 
1545  /* Make sure that arc "A" has the bigger radius */
1546  if ( radius_B > radius_A )
1547  {
1548  const POINT2D *tmp;
1549  tmp = B1; B1 = A1; A1 = tmp;
1550  tmp = B2; B2 = A2; A2 = tmp;
1551  tmp = B3; B3 = A3; A3 = tmp;
1552  P = CB; CB = CA; CA = P;
1553  d = radius_B; radius_B = radius_A; radius_A = d;
1554  }
1555 
1556  /* Circles touch at a point. Is that point within the arcs? */
1557  if ( d == (radius_A + radius_B) )
1558  {
1559  D.x = CA.x + (CB.x - CA.x) * radius_A / d;
1560  D.y = CA.y + (CB.y - CA.y) * radius_A / d;
1561 
1562  pt_in_arc_A = lw_pt_in_arc(&D, A1, A2, A3);
1563  pt_in_arc_B = lw_pt_in_arc(&D, B1, B2, B3);
1564 
1565  /* Arcs do touch at D, return it */
1566  if ( pt_in_arc_A && pt_in_arc_B )
1567  {
1568  dl->distance = 0.0;
1569  dl->p1 = D;
1570  dl->p2 = D;
1571  return LW_TRUE;
1572  }
1573  }
1574  /* Disjoint or contained circles don't intersect. Closest point may be on */
1575  /* the line joining CA to CB. */
1576  else if ( d > (radius_A + radius_B) /* Disjoint */ || d < (radius_A - radius_B) /* Contained */ )
1577  {
1578  POINT2D XA, XB; /* Points where the line from CA to CB cross their circle bounds */
1579 
1580  /* Calculate hypothetical nearest points, the places on the */
1581  /* two circles where the center-center line crosses. If both */
1582  /* arcs contain their hypothetical points, that's the crossing distance */
1583  XA.x = CA.x + (CB.x - CA.x) * radius_A / d;
1584  XA.y = CA.y + (CB.y - CA.y) * radius_A / d;
1585  XB.x = CB.x + (CA.x - CB.x) * radius_B / d;
1586  XB.y = CB.y + (CA.y - CB.y) * radius_B / d;
1587 
1588  pt_in_arc_A = lw_pt_in_arc(&XA, A1, A2, A3);
1589  pt_in_arc_B = lw_pt_in_arc(&XB, B1, B2, B3);
1590 
1591  /* If the nearest points are both within the arcs, that's our answer */
1592  /* the shortest distance is at the nearest points */
1593  if ( pt_in_arc_A && pt_in_arc_B )
1594  {
1595  return lw_dist2d_pt_pt(&XA, &XB, dl);
1596  }
1597  }
1598  /* Circles cross at two points, are either of those points in both arcs? */
1599  /* http://paulbourke.net/geometry/2circle/ */
1600  else if ( d < (radius_A + radius_B) )
1601  {
1602  POINT2D E, F; /* Points where circle(A) and circle(B) cross */
1603  /* Distance from CA to D */
1604  double a = (radius_A*radius_A - radius_B*radius_B + d*d) / (2*d);
1605  /* Distance from D to E or F */
1606  double h = sqrt(radius_A*radius_A - a*a);
1607 
1608  /* Location of D */
1609  D.x = CA.x + (CB.x - CA.x) * a / d;
1610  D.y = CA.y + (CB.y - CA.y) * a / d;
1611 
1612  /* Start from D and project h units perpendicular to CA-D to get E */
1613  E.x = D.x + (D.y - CA.y) * h / a;
1614  E.y = D.y + (D.x - CA.x) * h / a;
1615 
1616  /* Crossing point E contained in arcs? */
1617  pt_in_arc_A = lw_pt_in_arc(&E, A1, A2, A3);
1618  pt_in_arc_B = lw_pt_in_arc(&E, B1, B2, B3);
1619 
1620  if ( pt_in_arc_A && pt_in_arc_B )
1621  {
1622  dl->p1 = dl->p2 = E;
1623  dl->distance = 0.0;
1624  return LW_TRUE;
1625  }
1626 
1627  /* Start from D and project h units perpendicular to CA-D to get F */
1628  F.x = D.x - (D.y - CA.y) * h / a;
1629  F.y = D.y - (D.x - CA.x) * h / a;
1630 
1631  /* Crossing point F contained in arcs? */
1632  pt_in_arc_A = lw_pt_in_arc(&F, A1, A2, A3);
1633  pt_in_arc_B = lw_pt_in_arc(&F, B1, B2, B3);
1634 
1635  if ( pt_in_arc_A && pt_in_arc_B )
1636  {
1637  dl->p1 = dl->p2 = F;
1638  dl->distance = 0.0;
1639  return LW_TRUE;
1640  }
1641  }
1642  else
1643  {
1644  lwerror("lw_dist2d_arc_arc: arcs neither touch, intersect nor are disjoint! INCONCEIVABLE!");
1645  return LW_FALSE;
1646  }
1647 
1648  /* Closest point is in the arc A, but not in the arc B, so */
1649  /* one of the B end points must be the closest. */
1650  if (pt_in_arc_A && !pt_in_arc_B)
1651  {
1652  lw_dist2d_pt_arc(B1, A1, A2, A3, dl);
1653  lw_dist2d_pt_arc(B3, A1, A2, A3, dl);
1654  return LW_TRUE;
1655  }
1656  /* Closest point is in the arc B, but not in the arc A, so */
1657  /* one of the A end points must be the closest. */
1658  else if ( pt_in_arc_B && ! pt_in_arc_A )
1659  {
1660  lw_dist2d_pt_arc(A1, B1, B2, B3, dl);
1661  lw_dist2d_pt_arc(A3, B1, B2, B3, dl);
1662  return LW_TRUE;
1663  }
1664  /* Finally, one of the end-point to end-point combos is the closest. */
1665  else
1666  {
1667  lw_dist2d_pt_pt(A1, B1, dl);
1668  lw_dist2d_pt_pt(A1, B3, dl);
1669  lw_dist2d_pt_pt(A2, B1, dl);
1670  lw_dist2d_pt_pt(A2, B3, dl);
1671  return LW_TRUE;
1672  }
1673 
1674  return LW_TRUE;
1675 }
double lw_arc_center(const POINT2D *p1, const POINT2D *p2, const POINT2D *p3, POINT2D *result)
Determines the center of the circle defined by the three given points.
Definition: lwalgorithm.c:228
int lw_arc_is_pt(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3)
Returns true if arc A is actually a point (all vertices are the same) .
Definition: lwalgorithm.c:105
int lw_dist2d_arc_arc_concentric(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, double radius_A, const POINT2D *B1, const POINT2D *B2, const POINT2D *B3, double radius_B, const POINT2D *CENTER, DISTPTS *dl)
Definition: measures.c:1678
int mode
Definition: measures.h:51
POINT2D p1
Definition: measures.h:49
double x
Definition: liblwgeom.h:330
#define DIST_MIN
Definition: measures.h:41
#define LW_FALSE
Definition: liblwgeom.h:76
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:75
int lw_dist2d_pt_arc(const POINT2D *P, const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, DISTPTS *dl)
Definition: measures.c:1439
int lw_dist2d_seg_seg(const POINT2D *A, const POINT2D *B, const POINT2D *C, const POINT2D *D, DISTPTS *dl)
Finds the shortest distance between two segments.
Definition: measures.c:1820
POINT2D p2
Definition: measures.h:50
double y
Definition: liblwgeom.h:330
int lw_pt_in_arc(const POINT2D *P, const POINT2D *A1, const POINT2D *A2, const POINT2D *A3)
Returns true if P is on the same side of the plane partition defined by A1/A3 as A2 is...
Definition: lwalgorithm.c:85
double distance
Definition: measures.h:48
double distance2d_pt_pt(const POINT2D *p1, const POINT2D *p2)
Definition: measures.c:2313
int lw_dist2d_seg_arc(const POINT2D *A1, const POINT2D *A2, const POINT2D *B1, const POINT2D *B2, const POINT2D *B3, DISTPTS *dl)
Calculate the shortest distance between an arc and an edge.
Definition: measures.c:1292
#define FP_EQUALS(A, B)
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition: lwutil.c:190
int lw_dist2d_pt_pt(const POINT2D *thep1, const POINT2D *thep2, DISTPTS *dl)
Compares incoming points and stores the points closest to each other or most far away from each other...
Definition: measures.c:2281
Here is the call graph for this function:
Here is the caller graph for this function: