PostGIS  2.5.7dev-r@@SVN_REVISION@@

## ◆ lw_arc_center()

 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.

In the event the circle is complete, the midpoint of the segment defined by the first and second points is returned. If the points are collinear, as determined by equal slopes, then -1.0 is returned. If the interior point is coincident with either end point, they are taken as collinear. For non-collinear cases, arc radious is returned.

Definition at line 228 of file lwalgorithm.c.

229 {
230  POINT2D c;
231  double cx, cy, cr;
232  double dx21, dy21, dx31, dy31, h21, h31, d;
233
234  c.x = c.y = 0.0;
235
236  LWDEBUGF(2, "lw_arc_center called (%.16f,%.16f), (%.16f,%.16f), (%.16f,%.16f).", p1->x, p1->y, p2->x, p2->y, p3->x, p3->y);
237
238  /* Closed circle */
239  if (fabs(p1->x - p3->x) < EPSILON_SQLMM &&
240  fabs(p1->y - p3->y) < EPSILON_SQLMM)
241  {
242  cx = p1->x + (p2->x - p1->x) / 2.0;
243  cy = p1->y + (p2->y - p1->y) / 2.0;
244  c.x = cx;
245  c.y = cy;
246  *result = c;
247  cr = sqrt(pow(cx - p1->x, 2.0) + pow(cy - p1->y, 2.0));
248  return cr;
249  }
250
251  /* Using cartesian eguations from page https://en.wikipedia.org/wiki/Circumscribed_circle */
252  dx21 = p2->x - p1->x;
253  dy21 = p2->y - p1->y;
254  dx31 = p3->x - p1->x;
255  dy31 = p3->y - p1->y;
256
257  h21 = pow(dx21, 2.0) + pow(dy21, 2.0);
258  h31 = pow(dx31, 2.0) + pow(dy31, 2.0);
259
260  /* 2 * |Cross product|, d<0 means clockwise and d>0 counterclockwise sweeping angle */
261  d = 2 * (dx21 * dy31 - dx31 * dy21);
262
263  /* Check colinearity, |Cross product| = 0 */
264  if (fabs(d) < EPSILON_SQLMM)
265  return -1.0;
266
267  /* Calculate centroid coordinates and radius */
268  cx = p1->x + (h21 * dy31 - h31 * dy21) / d;
269  cy = p1->y - (h21 * dx31 - h31 * dx21) / d;
270  c.x = cx;
271  c.y = cy;
272  *result = c;
273  cr = sqrt(pow(cx - p1->x, 2) + pow(cy - p1->y, 2));
274
275  LWDEBUGF(2, "lw_arc_center center is (%.16f,%.16f)", result->x, result->y);
276
277  return cr;
278 }
#define EPSILON_SQLMM
Tolerance used to determine equality.
#define LWDEBUGF(level, msg,...)
Definition: lwgeom_log.h:88
double y
Definition: liblwgeom.h:331
double x
Definition: liblwgeom.h:331

References EPSILON_SQLMM, LWDEBUGF, POINT2D::x, and POINT2D::y.

Here is the caller graph for this function: