]> granicus.if.org Git - postgis/commitdiff
Change signatures for lwgeom distance sphere
authorPaul Ramsey <pramsey@cleverelephant.ca>
Tue, 6 Oct 2009 04:50:35 +0000 (04:50 +0000)
committerPaul Ramsey <pramsey@cleverelephant.ca>
Tue, 6 Oct 2009 04:50:35 +0000 (04:50 +0000)
git-svn-id: http://svn.osgeo.org/postgis/trunk@4599 b70326c6-7e19-0410-871a-916f4a2858ee

liblwgeom/cunit/cu_geodetic.c
liblwgeom/libgeom.h
liblwgeom/lwgeodetic.c
liblwgeom/lwgeodetic.h

index 825f4bfc7bb0a8f38d2cd248fcfc9facc278312f..95db1b848eac510c4d692b2da9277c8c96349701 100644 (file)
@@ -34,7 +34,6 @@ CU_pSuite register_geodetic_suite(void)
            (NULL == CU_add_test(pSuite, "test_clairaut()", test_clairaut))  || 
            (NULL == CU_add_test(pSuite, "test_edge_intersection()", test_edge_intersection))  ||
            (NULL == CU_add_test(pSuite, "test_edge_distance_to_point()", test_edge_distance_to_point)) ||
-           (NULL == CU_add_test(pSuite, "test_ptarray_point_in_ring_winding()", test_ptarray_point_in_ring_winding)) ||
            (NULL == CU_add_test(pSuite, "test_edge_distance_to_edge()", test_edge_distance_to_edge)) || 
            (NULL == CU_add_test(pSuite, "test_lwgeom_distance_sphere()", test_lwgeom_distance_sphere)) ||
            (NULL == CU_add_test(pSuite, "test_ptarray_point_in_ring()", test_ptarray_point_in_ring)) 
@@ -487,7 +486,9 @@ void test_lwgeom_distance_sphere(void)
        /* Line/line distance, 1 degree apart */
        lwg1 = lwgeom_from_ewkt("LINESTRING(-30 10, -20 5, -10 3, 0 1)", PARSER_CHECK_NONE);
        lwg2 = lwgeom_from_ewkt("LINESTRING(-10 -5, -5 0, 5 0, 10 -5)", PARSER_CHECK_NONE);
-       d = lwgeom_distance_sphere(lwg1, lwg2, NULL, NULL, 0.0);        
+       lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
+       lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, M_PI / 180.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -495,7 +496,9 @@ void test_lwgeom_distance_sphere(void)
        /* Line/line distance, crossing, 0.0 apart */
        lwg1 = lwgeom_from_ewkt("LINESTRING(-30 10, -20 5, -10 3, 0 1)", PARSER_CHECK_NONE);
        lwg2 = lwgeom_from_ewkt("LINESTRING(-10 -5, -5 20, 5 0, 10 -5)", PARSER_CHECK_NONE);
-       d = lwgeom_distance_sphere(lwg1, lwg2, NULL, NULL, 0.0);        
+       lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
+       lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, 0.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -503,14 +506,18 @@ void test_lwgeom_distance_sphere(void)
        /* Line/point distance, 1 degree apart */
        lwg1 = lwgeom_from_ewkt("POINT(-4 1)", PARSER_CHECK_NONE);
        lwg2 = lwgeom_from_ewkt("LINESTRING(-10 -5, -5 0, 5 0, 10 -5)", PARSER_CHECK_NONE);
-       d = lwgeom_distance_sphere(lwg1, lwg2, NULL, NULL, 0.0);        
+       lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
+       lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, M_PI / 180.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
 
        lwg1 = lwgeom_from_ewkt("POINT(-4 1)", PARSER_CHECK_NONE);
        lwg2 = lwgeom_from_ewkt("POINT(-4 -1)", PARSER_CHECK_NONE);
-       d = lwgeom_distance_sphere(lwg1, lwg2, NULL, NULL, 0.0);        
+       lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
+       lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, M_PI / 90.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -520,7 +527,7 @@ void test_lwgeom_distance_sphere(void)
        lwg2 = lwgeom_from_ewkt("POINT(-1 -1)", PARSER_CHECK_NONE);
        lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
        lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
-       d = lwgeom_distance_sphere(lwg1, lwg2, &gbox1, &gbox2, 0.0);    
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, 0.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -530,7 +537,7 @@ void test_lwgeom_distance_sphere(void)
        lwg2 = lwgeom_from_ewkt("POINT(-1 -1)", PARSER_CHECK_NONE);
        lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
        lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
-       d = lwgeom_distance_sphere(lwg1, lwg2, &gbox1, &gbox2, 0.0);    
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, M_PI / 180.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -540,7 +547,7 @@ void test_lwgeom_distance_sphere(void)
        lwg2 = lwgeom_from_ewkt("POINT(2 2)", PARSER_CHECK_NONE);
        lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
        lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
-       d = lwgeom_distance_sphere(lwg1, lwg2, &gbox1, &gbox2, 0.0);    
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);      
        CU_ASSERT_DOUBLE_EQUAL(d, 0.0, 0.00001);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
@@ -550,7 +557,7 @@ void test_lwgeom_distance_sphere(void)
        lwg2 = lwgeom_from_ewkt("0106000020E61000000100000001030000000100000007000000280EC3FB8CCA5EC0A5CDC747233C45402787C8F58CCA5EC0659EA2761E3C45400CED58DF8FCA5EC0C37FAE6E1E3C4540AE97B8E08FCA5EC00346F58B1F3C4540250359FD8ECA5EC05460628E1F3C45403738F4018FCA5EC05DC84042233C4540280EC3FB8CCA5EC0A5CDC747233C4540", PARSER_CHECK_NONE);
        lwgeom_calculate_gbox_geodetic(lwg1, &gbox1);
        lwgeom_calculate_gbox_geodetic(lwg2, &gbox2);
-       d = lwgeom_distance_sphere(lwg1, lwg2, &gbox1, &gbox2, 0.0);
+       d = lwgeom_distance_sphere(lwg1, lwg2, gbox1, gbox2, 0.0);
        CU_ASSERT_DOUBLE_EQUAL(d * 6371009.0, 23630.8003, 0.1);
        lwgeom_free(lwg1);
        lwgeom_free(lwg2);
index ff1602f8137df04ad068057f12c58ffae98ef9e1..ada2c4a35da09dd5002f2ea6a45f71fdef7ffc36 100644 (file)
@@ -374,7 +374,7 @@ extern int lwgeom_calculate_gbox(const LWGEOM *lwgeom, GBOX *gbox);
 * Calculate the geodetic distance from lwgeom1 to lwgeom2 on the unit sphere. 
 * Pass in a tolerance in radians.
 */
-extern double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBOX *gbox2, double tolerance);
+extern double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX gbox1, GBOX gbox2, double tolerance);
 
 /**
 * New function to read doubles directly from the double* coordinate array
index 599f9c4ee09a75a05e04ed241f0997c575d912f0..214f63345dfb8101a1f5006d275ac0376953c23f 100644 (file)
@@ -1441,7 +1441,7 @@ static double ptarray_distance_sphere(POINTARRAY *pa1, POINTARRAY *pa2, double t
 * longitude and latitude. Return immediately when the calulated distance drops
 * below the tolerance (useful for dwithin calculations).
 */
-double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBOX *gbox2, double tolerance)
+double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX gbox1, GBOX gbox2, double tolerance)
 {
        int type1, type2;
        
@@ -1478,12 +1478,6 @@ double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBO
                return ptarray_distance_sphere(pa1, pa2, tolerance);
        }
        
-       if( ! gbox1 || ! gbox2 ) 
-       {
-               lwerror("gboxes are required to calculate distances from spherical lwgeoms");
-               return -1.0;
-       }
-       
        /* Point/Polygon cases, if point-in-poly, return zero, else return distance. */
        if( ( type1 == POLYGONTYPE && type2 == POINTTYPE ) || 
            ( type2 == POLYGONTYPE && type1 == POINTTYPE ) )
@@ -1491,7 +1485,7 @@ double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBO
                POINT2D p;
                LWPOLY *lwpoly;
                LWPOINT *lwpt;
-               GBOX *gbox;
+               GBOX gbox;
                double distance = MAXFLOAT;
                int i;
                
@@ -1532,7 +1526,7 @@ double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBO
                POINT2D p;
                LWPOLY *lwpoly;
                LWLINE *lwline;
-               GBOX *gbox;
+               GBOX gbox;
                double distance = MAXFLOAT;
                int i;
                
@@ -1655,19 +1649,13 @@ double lwgeom_distance_sphere(LWGEOM *lwgeom1, LWGEOM *lwgeom2, GBOX *gbox1, GBO
 * a guaranteed outside point (lon/lat decimal degrees) (calculate with gbox_pt_outside()) 
 * return LW_TRUE if point is inside or on edge of polygon.
 */
-int lwpoly_covers_point2d(const LWPOLY *poly, GBOX *gbox, POINT2D pt_to_test)
+int lwpoly_covers_point2d(const LWPOLY *poly, GBOX gbox, POINT2D pt_to_test)
 {
        int i;
        int in_hole_count = 0;
-       POINT3D p;
-       GEOGRAPHIC_POINT g;
+       POINT3D p, q;
+       GEOGRAPHIC_POINT g, gpt_to_test;
        POINT2D pt_outside;
-
-       if( ! gbox )
-       {
-               lwerror("gbox is required to calculate spherical point-in-poly with lwgeom");
-               return LW_FALSE;
-       }
        
        /* Nulls and empties don't contain anything! */
        if( ! poly || lwgeom_is_empty((LWGEOM*)poly) )
@@ -1676,8 +1664,14 @@ int lwpoly_covers_point2d(const LWPOLY *poly, GBOX *gbox, POINT2D pt_to_test)
                return LW_FALSE;
        }
 
+       /* Point not in box? Done! */
+       geographic_point_init(pt_to_test.x, pt_to_test.y, &gpt_to_test);
+       geog2cart(gpt_to_test, &q);
+       if( ! gbox_contains_point3d(gbox, q) )
+               return LW_FALSE;
+       
        /* Calculate our outside point from the gbox */
-       gbox_pt_outside(*gbox, &p);
+       gbox_pt_outside(gbox, &p);
        cart2geog(p, &g);
        pt_outside.x = rad2deg(g.lon);
        pt_outside.y = rad2deg(g.lat);
index 8577a05586105e62d7393a07487c657f04cd76a4..153c48687a7677d1492f09c849d91fbc67043c23 100644 (file)
@@ -73,5 +73,5 @@ void point_deg2rad(GEOGRAPHIC_POINT *p);
 void point_rad2deg(GEOGRAPHIC_POINT *p);
 void geographic_point_init(double lon, double lat, GEOGRAPHIC_POINT *g);
 int ptarray_point_in_ring_winding(POINTARRAY *pa, POINT2D pt_to_test);
-int lwpoly_covers_point2d(const LWPOLY *poly, GBOX *gbox, POINT2D pt_to_test);
+int lwpoly_covers_point2d(const LWPOLY *poly, GBOX gbox, POINT2D pt_to_test);
 int ptarray_point_in_ring(POINTARRAY *pa, POINT2D pt_outside, POINT2D pt_to_test);