diff --git a/liblwgeom/cunit/cu_ptarray.c b/liblwgeom/cunit/cu_ptarray.c index b1428dd95..3e562a02f 100644 --- a/liblwgeom/cunit/cu_ptarray.c +++ b/liblwgeom/cunit/cu_ptarray.c @@ -472,7 +472,7 @@ static void test_ptarray_contains_point_arc() POINTARRAY *pa; POINT2D pt; int rv; - + /* Collection of semi-circles surrounding unit square */ lwline = lwgeom_as_lwline(lwgeom_from_text("LINESTRING(-1 -1, -2 0, -1 1, 0 2, 1 1, 2 0, 1 -1, 0 -2, -1 -1)")); pa = lwline->points; @@ -509,18 +509,25 @@ static void test_ptarray_contains_point_arc() /* Two-edge ring made up of semi-circles (really, a circle) */ lwfree(lwline); + lwline = lwgeom_as_lwline(lwgeom_from_text("LINESTRING(-1 0, 0 1, 1 0, 0 -1, -1 0)")); pa = lwline->points; - + /* Point outside */ pt.x = -1.5; pt.y = 1.5; rv = ptarray_contains_point_arc(pa, &pt); CU_ASSERT_EQUAL(rv, -1); - /* Point inside */ - pt.x = -0.2; - pt.y = 0.2; + /* Point inside at middle */ + pt.x = 0; + pt.y = 0; + rv = ptarray_contains_point_arc(pa, &pt); + CU_ASSERT_EQUAL(rv, 1); + + /* Point inside offset from middle */ + pt.x = 0.01; + pt.y = 0.01; rv = ptarray_contains_point_arc(pa, &pt); CU_ASSERT_EQUAL(rv, 1); @@ -559,7 +566,8 @@ static void test_ptarray_contains_point_arc() pa = lwline->points; cu_error_msg_reset(); rv = ptarray_contains_point_arc(pa, &pt); - CU_ASSERT_STRING_EQUAL("ptarray_contains_point_arc called on unclosed ring", cu_error_msg); + //printf("%s\n", cu_error_msg); + CU_ASSERT_STRING_EQUAL("ptarray_contains_point_arc called with even number of points", cu_error_msg); /* Unclosed ring */ lwfree(lwline); @@ -571,6 +579,7 @@ static void test_ptarray_contains_point_arc() lwline_free(lwline); } + /* ** Used by the test harness to register the tests in this file. */ diff --git a/liblwgeom/liblwgeom_internal.h b/liblwgeom/liblwgeom_internal.h index 43404a090..240580904 100644 --- a/liblwgeom/liblwgeom_internal.h +++ b/liblwgeom/liblwgeom_internal.h @@ -353,6 +353,7 @@ int lwtin_is_closed(const LWTIN *tin); * Returns -1 for left and 1 for right and 0 for co-linearity */ int lw_segment_side(const POINT2D *p1, const POINT2D *p2, const POINT2D *q); +int lw_arc_side(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, const POINT2D *Q); int lw_arc_calculate_gbox_cartesian_2d(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, GBOX *gbox); double lw_arc_center(const POINT2D *p1, const POINT2D *p2, const POINT2D *p3, POINT2D *result); int lw_pt_in_seg(const POINT2D *P, const POINT2D *A1, const POINT2D *A2); @@ -360,7 +361,6 @@ int lw_pt_in_arc(const POINT2D *P, const POINT2D *A1, const POINT2D *A2, const P int lw_arc_is_pt(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3); double lw_seg_length(const POINT2D *A1, const POINT2D *A2); double lw_arc_length(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3); -double lw_arc_side(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, const POINT2D *Q); int pt_in_ring_2d(const POINT2D *p, const POINTARRAY *ring); int ptarray_contains_point(const POINTARRAY *pa, const POINT2D *pt); int ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt); diff --git a/liblwgeom/lwalgorithm.c b/liblwgeom/lwalgorithm.c index 4dff6f867..2768cb073 100644 --- a/liblwgeom/lwalgorithm.c +++ b/liblwgeom/lwalgorithm.c @@ -174,7 +174,7 @@ lw_arc_length(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3) return circumference_A * (angle / (2*M_PI)); } -double lw_arc_side(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, const POINT2D *Q) +int lw_arc_side(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, const POINT2D *Q) { POINT2D C; double radius_A; @@ -197,6 +197,12 @@ double lw_arc_side(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, cons return 0; } + /* Q on A1-A3 line, so its on opposite side to A2 */ + if ( side_Q == 0 ) + { + return -1 * side_A2; + } + /* * Q is inside the arc boundary, so it's not on the side we * might think from examining only the end points diff --git a/liblwgeom/measures.c b/liblwgeom/measures.c index dbee74e90..382d6863a 100644 --- a/liblwgeom/measures.c +++ b/liblwgeom/measures.c @@ -665,7 +665,6 @@ lw_dist2d_poly_poly(LWPOLY *poly1, LWPOLY *poly2, DISTPTS *dl) } /** - * search all the segments of pointarray to see which one is closest to p1 * Returns minimum distance between point and pointarray */ @@ -695,6 +694,47 @@ lw_dist2d_pt_ptarray(POINT2D *p, POINTARRAY *pa,DISTPTS *dl) return LW_TRUE; } +/** + * search all the segments of pointarray to see which one is closest to p1 + * Returns minimum distance between point and pointarray + */ +#if 0 +int +lw_dist2d_pt_ptarrayarc(const POINT2D *p, const POINTARRAY *pa, DISTPTS *dl) +{ + int t; + const POINT2D *A1; + const POINT2D *A2; + const POINT2D *A3; + int twist = dl->twisted; + + LWDEBUG(2, "lw_dist2d_pt_ptarrayarc is called"); + + A1 = getPoint2d_cp(pa, 0); + A2 = getPoint2d_cp(pa, 1); + + if ( !lw_dist2d_pt_pt(p, A1, dl) ) + return LW_FALSE; + + for ( t=2; tnpoints; t++ ) + { + dl->twisted=twist; + A3 = getPoint2d_cp(pa, t); + + if ( lw_dist2d_pt_arc(p, A1, A2, A3, dl) == LW_FALSE ) + return LW_FALSE; + + if ( dl->distance<=dl->tolerance && dl->mode == DIST_MIN ) + return LW_TRUE; /*just a check if the answer is already given*/ + + A1 = A2; + A2 = A3; + } + + return LW_TRUE; +} +#endif + /** * Brute force. * Test line-ring distance against each ring. diff --git a/liblwgeom/ptarray.c b/liblwgeom/ptarray.c index 0b7642668..86c890403 100644 --- a/liblwgeom/ptarray.c +++ b/liblwgeom/ptarray.c @@ -793,11 +793,17 @@ ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt) /* Check for not an arc ring (always have odd # of points) */ if ( (pa->npoints % 2) == 0 ) + { lwerror("ptarray_contains_point_arc called with even number of points"); + return -1; + } /* Check for not an arc ring (always have >= 3 points) */ if ( pa->npoints < 3 ) + { lwerror("ptarray_contains_point_arc called too-short pointarray"); + return -1; + } /* Check for unclosed case */ seg1 = getPoint2d_cp(pa, 0); @@ -805,6 +811,7 @@ ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt) if ( ! p2d_same(seg1, seg3) ) { lwerror("ptarray_contains_point_arc called on unclosed ring"); + return -1; } /* OK, it's closed. Is it just one circle? */ else if ( pa->npoints == 3 ) @@ -833,10 +840,11 @@ ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt) } /* Start on the ring */ - seg2 = getPoint2d_cp(pa, 1); - for ( i=2; i < pa->npoints; i++ ) + seg1 = getPoint2d_cp(pa, 0); + for ( i=1; i < pa->npoints; i += 2 ) { - seg3 = getPoint2d_cp(pa, i); + seg2 = getPoint2d_cp(pa, i); + seg3 = getPoint2d_cp(pa, i+1); /* Catch an easy boundary case */ if( p2d_same(seg3, pt) ) @@ -845,18 +853,15 @@ ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt) /* Skip arcs that have no size */ if ( lw_arc_is_pt(seg1, seg2, seg3) ) { - seg1 = seg2; - seg2 = seg3; + seg1 = seg3; continue; } - /* Only test segments in our vertical range */ lw_arc_calculate_gbox_cartesian_2d(seg1, seg2, seg3, &gbox); if ( pt->y > gbox.ymax || pt->y < gbox.ymin ) { - seg1 = seg2; - seg2 = seg3; + seg1 = seg3; continue; } @@ -880,8 +885,30 @@ ptarray_contains_point_arc(const POINTARRAY *pa, const POINT2D *pt) wn--; } - seg1 = seg2; - seg2 = seg3; + /* Inside the arc! */ + if ( pt->x <= gbox.xmax && pt->x >= gbox.xmin ) + { + POINT2D C; + double radius = lw_arc_center(seg1, seg2, seg3, &C); + double d = distance2d_pt_pt(pt, &C); + + /* On the boundary! */ + if ( d == radius ) + return 0; + + /* Within the arc! */ + if ( d < radius ) + { + /* Left side, increment winding number */ + if ( side < 0 ) + wn++; + /* Right side, decrement winding number */ + if ( side > 0 ) + wn--; + } + } + + seg1 = seg3; } /* Outside */