From: Paul Ramsey Date: Sat, 4 Feb 2012 01:30:14 +0000 (+0000) Subject: Geography should support SRIDs other than 4326 (#1538) X-Git-Tag: 2.0.0alpha4~75 X-Git-Url: https://granicus.if.org/sourcecode?a=commitdiff_plain;h=c1f7407a09397c45338db29f105628cda655449b;p=postgis Geography should support SRIDs other than 4326 (#1538) git-svn-id: http://svn.osgeo.org/postgis/trunk@9029 b70326c6-7e19-0410-871a-916f4a2858ee --- diff --git a/liblwgeom/liblwgeom.h.in b/liblwgeom/liblwgeom.h.in index 067d384b3..a5dfad141 100644 --- a/liblwgeom/liblwgeom.h.in +++ b/liblwgeom/liblwgeom.h.in @@ -19,6 +19,7 @@ #include #include #include +#include "projects.h" #include "proj_api.h" /** @@ -133,6 +134,12 @@ #define SRID_UNKNOWN 0 #define SRID_IS_UNKNOWN(x) ((int)x<=0) +/* +** EPSG WGS84 geographics, OGC standard default SRS, better be in +** the SPATIAL_REF_SYS table! +*/ +#define SRID_DEFAULT 4326 + /** * Return a valid SRID from an arbitrary integer * Raises a notice if what comes out is different from diff --git a/liblwgeom/liblwgeom_internal.h b/liblwgeom/liblwgeom_internal.h index 26f533641..481646301 100644 --- a/liblwgeom/liblwgeom_internal.h +++ b/liblwgeom/liblwgeom_internal.h @@ -34,7 +34,7 @@ /** * PI */ -#define PI 3.1415926535897932384626433832795 +#define LW_PI 3.1415926535897932384626433832795 /** * Floating point comparitors. diff --git a/liblwgeom/lwgeom_transform.c b/liblwgeom/lwgeom_transform.c index b68676766..76884c466 100644 --- a/liblwgeom/lwgeom_transform.c +++ b/liblwgeom/lwgeom_transform.c @@ -10,7 +10,6 @@ * **********************************************************************/ -#include "proj_api.h" #include "liblwgeom.h" #include "lwgeom_log.h" #include diff --git a/libpgcommon/lwgeom_transform.c b/libpgcommon/lwgeom_transform.c index 35b22afe7..950061456 100644 --- a/libpgcommon/lwgeom_transform.c +++ b/libpgcommon/lwgeom_transform.c @@ -10,30 +10,32 @@ * **********************************************************************/ + +/* PostgreSQL headers */ #include "postgres.h" #include "fmgr.h" #include "miscadmin.h" +#include "utils/memutils.h" +#include "executor/spi.h" +#include "access/hash.h" +#include "utils/hsearch.h" +/* PostGIS headers */ #include "liblwgeom.h" #include "lwgeom_pg.h" +#include "lwgeom_transform.h" +/* C headers */ #include #include #include #include -#include "proj_api.h" -#include "utils/memutils.h" -#include "executor/spi.h" -#include "access/hash.h" -#include "utils/hsearch.h" -#include "lwgeom_transform.h" - +/* Expose an internal Proj function */ int pj_transform_nodatum(projPJ srcdefn, projPJ dstdefn, long point_count, int point_offset, double *x, double *y, double *z ); - /* PROJ 4 lookup transaction cache methods */ #define PROJ4_CACHE_ITEMS 8 @@ -95,7 +97,7 @@ static projPJ GetPJHashEntry(MemoryContext mcxt); static void DeletePJHashEntry(MemoryContext mcxt); /* Internal Cache API */ -static PROJ4PortalCache *GetPROJ4SRSCache(FunctionCallInfoData *fcinfo) ; +static PROJ4PortalCache *GetPROJ4SRSCache(FunctionCallInfo fcinfo) ; static bool IsInPROJ4SRSCache(PROJ4PortalCache *PROJ4Cache, int srid); static projPJ GetProjectionFromPROJ4SRSCache(PROJ4PortalCache *PROJ4Cache, int srid); static void AddToPROJ4SRSCache(PROJ4PortalCache *PROJ4Cache, int srid, int other_srid); @@ -626,13 +628,11 @@ void SetPROJ4LibPath(void) } } - - -Proj4Cache GetPROJ4Cache(FunctionCallInfoData *fcinfo) { +Proj4Cache GetPROJ4Cache(FunctionCallInfo fcinfo) { return (Proj4Cache)GetPROJ4SRSCache(fcinfo) ; } -static PROJ4PortalCache *GetPROJ4SRSCache(FunctionCallInfoData *fcinfo) +static PROJ4PortalCache *GetPROJ4SRSCache(FunctionCallInfo fcinfo) { PROJ4PortalCache *PROJ4Cache ; @@ -677,7 +677,77 @@ static PROJ4PortalCache *GetPROJ4SRSCache(FunctionCallInfoData *fcinfo) } +int +GetProjectionsUsingFCInfo(FunctionCallInfo fcinfo, int srid1, int srid2, projPJ *pj1, projPJ *pj2) +{ + Proj4Cache *proj_cache = NULL; + /* Set the search path if we haven't already */ + SetPROJ4LibPath(); + /* get or initialize the cache for this round */ + proj_cache = GetPROJ4Cache(fcinfo); + if ( !proj_cache ) + return LW_FAILURE; + + /* Add the output srid to the cache if it's not already there */ + if (!IsInPROJ4Cache(proj_cache, srid1)) + AddToPROJ4Cache(proj_cache, srid1, srid2); + + /* Add the input srid to the cache if it's not already there */ + if (!IsInPROJ4Cache(proj_cache, srid2)) + AddToPROJ4Cache(proj_cache, srid2, srid1); + + /* Get the projections */ + *pj1 = GetProjectionFromPROJ4Cache(proj_cache, srid1); + *pj2 = GetProjectionFromPROJ4Cache(proj_cache, srid2); + + return LW_SUCCESS; +} + +int +spheroid_init_from_srid(FunctionCallInfo fcinfo, int srid, SPHEROID *s) +{ + projPJ pj1; + projPJ pj2; + PJ *p; + + if ( GetProjectionsUsingFCInfo(fcinfo, srid, srid, &pj1, &pj2) == LW_FAILURE) + return LW_FAILURE; + + if ( ! pj_is_latlong(pj1) ) + return LW_FAILURE; + + /* Get access to the proj internals */ + p = (PJ*)pj1; + + /* Initialize */ + s->a = p->a; + s->e_sq = p->es; + s->b = s->a * sqrt(p->one_es); + s->f = (s->a - s->b) / s->a; + s->radius = (2.0 * s->a + s->b ) / 3.0; + + return LW_SUCCESS; +} + +void srid_is_latlong(FunctionCallInfo fcinfo, int srid) +{ + projPJ pj1; + projPJ pj2; + + if ( srid == SRID_DEFAULT || srid == SRID_UNKNOWN ) + return; + + if ( GetProjectionsUsingFCInfo(fcinfo, srid, srid, &pj1, &pj2) == LW_FAILURE) + return; + + if ( pj_is_latlong(pj1) ) + return; + + ereport(ERROR, ( + errcode(ERRCODE_INVALID_PARAMETER_VALUE), + errmsg("Only lon/lat coordinate systems are supported in geography."))); +} diff --git a/libpgcommon/lwgeom_transform.h b/libpgcommon/lwgeom_transform.h index 0b67623cf..235151f57 100644 --- a/libpgcommon/lwgeom_transform.h +++ b/libpgcommon/lwgeom_transform.h @@ -13,7 +13,8 @@ #include "postgres.h" #include "liblwgeom.h" #include "lwgeom_pg.h" -#include "proj_api.h" + + char* GetProj4StringSPI(int srid); void SetPROJ4LibPath(void) ; @@ -25,12 +26,14 @@ void SetPROJ4LibPath(void) ; typedef void *Proj4Cache ; void SetPROJ4LibPath(void); -Proj4Cache GetPROJ4Cache(FunctionCallInfoData *fcinfo) ; +Proj4Cache GetPROJ4Cache(FunctionCallInfo fcinfo) ; bool IsInPROJ4Cache(Proj4Cache cache, int srid) ; void AddToPROJ4Cache(Proj4Cache cache, int srid, int other_srid); void DeleteFromPROJ4Cache(Proj4Cache cache, int srid) ; projPJ GetProjectionFromPROJ4Cache(Proj4Cache cache, int srid); - +int GetProjectionsUsingFCInfo(FunctionCallInfo fcinfo, int srid1, int srid2, projPJ *pj1, projPJ *pj2); +int spheroid_init_from_srid(FunctionCallInfo fcinfo, int srid, SPHEROID *s); +void srid_is_latlong(FunctionCallInfo fcinfo, int srid); /** * Builtin SRID values diff --git a/postgis/geography.h b/postgis/geography.h index ceea39ea3..cf07dd4fa 100644 --- a/postgis/geography.h +++ b/postgis/geography.h @@ -23,12 +23,6 @@ #define WGS84_MINOR_AXIS (WGS84_MAJOR_AXIS - WGS84_MAJOR_AXIS / WGS84_INVERSE_FLATTENING) #define WGS84_RADIUS ((2.0 * WGS84_MAJOR_AXIS + WGS84_MINOR_AXIS ) / 3.0) -/* -** EPSG WGS84 geographics, OGC standard default SRS, better be in -** the SPATIAL_REF_SYS table! -*/ -#define SRID_DEFAULT 4326 - /********************************************************************** ** Useful functions for all GSERIALIZED handlers. ** TODO: Move to common.h in pgcommon diff --git a/postgis/geography_inout.c b/postgis/geography_inout.c index d121a4b21..4583a5d68 100644 --- a/postgis/geography_inout.c +++ b/postgis/geography_inout.c @@ -29,6 +29,7 @@ #include "lwgeom_pg.h" /* For debugging macros. */ #include "geography.h" /* For utility functions. */ #include "lwgeom_export.h" /* For export functions. */ +#include "lwgeom_transform.h" Datum geography_in(PG_FUNCTION_ARGS); Datum geography_out(PG_FUNCTION_ARGS); @@ -161,7 +162,10 @@ Datum geography_in(PG_FUNCTION_ARGS) lwgeom = lwg_parser_result.geom; } + /* Error on any SRID != default */ + srid_is_latlong(fcinfo, lwgeom->srid); + /* Convert to gserialized */ g_ser = gserialized_geography_from_lwgeom(lwgeom, geog_typmod); /* Clean up temporary object */ @@ -500,14 +504,25 @@ Datum geography_as_geojson(PG_FUNCTION_ARGS) PG_FUNCTION_INFO_V1(geography_from_text); Datum geography_from_text(PG_FUNCTION_ARGS) { + LWGEOM_PARSER_RESULT lwg_parser_result; + GSERIALIZED *g_ser = NULL; text *wkt_text = PG_GETARG_TEXT_P(0); + /* Extract the cstring from the varlena */ char *wkt = text2cstring(wkt_text); - /* Pass the cstring to the input parser, and magic occurs! */ - Datum rv = DirectFunctionCall3(geography_in, PointerGetDatum(wkt), Int32GetDatum(0), Int32GetDatum(-1)); - /* Clean up and return */ + + /* Pass the cstring to the input parser, and magic occurs! */ + if ( lwgeom_parse_wkt(&lwg_parser_result, wkt, LW_PARSER_CHECK_ALL) == LW_FAILURE ) + PG_PARSER_ERROR(lwg_parser_result); + + /* Clean up string */ pfree(wkt); - PG_RETURN_DATUM(rv); + g_ser = gserialized_geography_from_lwgeom(lwg_parser_result.geom, 0); + + /* Clean up temporary object */ + lwgeom_free(lwg_parser_result.geom); + + PG_RETURN_POINTER(g_ser); } /* @@ -516,27 +531,21 @@ Datum geography_from_text(PG_FUNCTION_ARGS) PG_FUNCTION_INFO_V1(geography_from_binary); Datum geography_from_binary(PG_FUNCTION_ARGS) { - char *wkb_cstring = NULL; - text *wkb_hex = NULL; char *wkb_bytea = (char*)PG_DETOAST_DATUM(PG_GETARG_DATUM(0)); - char *hexarg = palloc(4 + VARHDRSZ); - - /* Create our second argument to binary_encode */ - SET_VARSIZE(hexarg, 4 + VARHDRSZ); - memcpy(VARDATA(hexarg), "hex", 4); - - /* Convert the bytea into a hex representation cstring */ - wkb_hex = (text*)DatumGetPointer(DirectFunctionCall2(binary_encode, PointerGetDatum(wkb_bytea), PointerGetDatum(hexarg))); - wkb_cstring = text2cstring(wkb_hex); - pfree(hexarg); - - /* Pass the cstring to the input parser, and magic occurs! */ - PG_RETURN_DATUM(DirectFunctionCall3(geography_in, PointerGetDatum(wkb_cstring), Int32GetDatum(0), Int32GetDatum(-1))); + GSERIALIZED *gser = NULL; + size_t wkb_size = VARSIZE(wkb_bytea); + uint8_t *wkb = (uint8_t*)VARDATA(wkb_bytea); + LWGEOM *lwgeom = lwgeom_from_wkb(wkb, wkb_size, LW_PARSER_CHECK_NONE); + + if ( ! lwgeom ) + lwerror("Unable to parse WKB"); + + gser = gserialized_geography_from_lwgeom(lwgeom, 0); + lwgeom_free(lwgeom); + PG_RETURN_POINTER(gser); } - - PG_FUNCTION_INFO_V1(geography_from_geometry); Datum geography_from_geometry(PG_FUNCTION_ARGS) { @@ -555,12 +564,7 @@ Datum geography_from_geometry(PG_FUNCTION_ARGS) } /* Error on any SRID != default */ - if ( lwgeom->srid != SRID_DEFAULT ) - { - ereport(ERROR, ( - errcode(ERRCODE_INVALID_PARAMETER_VALUE), - errmsg("Only SRID %d is currently supported in geography.", SRID_DEFAULT))); - } + srid_is_latlong(fcinfo, lwgeom->srid); /* Check if the geography has valid coordinate range. */ if ( lwgeom_check_geodetic(lwgeom) == LW_FALSE ) @@ -632,6 +636,9 @@ Datum geography_recv(PG_FUNCTION_ARGS) lwgeom = lwgeom_from_wkb((uint8_t*)buf->data, buf->len, LW_PARSER_CHECK_ALL); + /* Error on any SRID != default */ + srid_is_latlong(fcinfo, lwgeom->srid); + g_ser = gserialized_geography_from_lwgeom(lwgeom, geog_typmod); /* Clean up temporary object */ diff --git a/postgis/geography_measurement.c b/postgis/geography_measurement.c index a81b4a07a..e31a4d35f 100644 --- a/postgis/geography_measurement.c +++ b/postgis/geography_measurement.c @@ -24,7 +24,7 @@ #include "liblwgeom_internal.h" /* For FP comparators. */ #include "lwgeom_pg.h" /* For debugging macros. */ #include "geography.h" /* For utility functions. */ -#include "lwgeom_transform.h" +#include "lwgeom_transform.h" /* For SRID functions */ Datum geography_distance(PG_FUNCTION_ARGS); Datum geography_dwithin(PG_FUNCTION_ARGS); @@ -65,8 +65,8 @@ Datum geography_distance(PG_FUNCTION_ARGS) use_spheroid = PG_GETARG_BOOL(3); /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); - + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s); + /* Set to sphere if requested */ if ( ! use_spheroid ) s.a = s.b = s.radius; @@ -129,7 +129,7 @@ Datum geography_dwithin(PG_FUNCTION_ARGS) use_spheroid = PG_GETARG_BOOL(3); /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s); /* Set to sphere if requested */ if ( ! use_spheroid ) @@ -222,7 +222,7 @@ Datum geography_area(PG_FUNCTION_ARGS) use_spheroid = PG_GETARG_BOOL(1); /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s); lwgeom = lwgeom_from_gserialized(g); @@ -311,7 +311,7 @@ Datum geography_perimeter(PG_FUNCTION_ARGS) use_spheroid = PG_GETARG_BOOL(1); /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s); /* User requests spherical calculation, turn our spheroid into a sphere */ if ( ! use_spheroid ) @@ -362,7 +362,7 @@ Datum geography_length(PG_FUNCTION_ARGS) use_spheroid = PG_GETARG_BOOL(1); /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s); /* User requests spherical calculation, turn our spheroid into a sphere */ if ( ! use_spheroid ) @@ -647,7 +647,7 @@ Datum geography_project(PG_FUNCTION_ARGS) azimuth = PG_GETARG_FLOAT8(2); /* Azimuth in Radians */ /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g), &s); /* Handle the zero distance case */ if( FP_EQUALS(distance, 0.0) ) @@ -717,7 +717,7 @@ Datum geography_azimuth(PG_FUNCTION_ARGS) } /* Initialize spheroid */ - spheroid_init(&s, WGS84_MAJOR_AXIS, WGS84_MINOR_AXIS); + spheroid_init_from_srid(fcinfo, gserialized_get_srid(g1), &s); /* Calculate the direction */ azimuth = lwgeom_azumith_spheroid(lwgeom_as_lwpoint(lwgeom1), lwgeom_as_lwpoint(lwgeom2), &s); diff --git a/postgis/lwgeom_transform.c b/postgis/lwgeom_transform.c index 91d6369cf..a006fdad0 100644 --- a/postgis/lwgeom_transform.c +++ b/postgis/lwgeom_transform.c @@ -12,7 +12,6 @@ #include "postgres.h" #include "fmgr.h" -#include "proj_api.h" #include "liblwgeom.h" #include "lwgeom_transform.h" @@ -22,33 +21,6 @@ Datum transform(PG_FUNCTION_ARGS); Datum transform_geom(PG_FUNCTION_ARGS); Datum postgis_proj_version(PG_FUNCTION_ARGS); -static int -GetProjectionsUsingFCInfo(FunctionCallInfo fcinfo, int srid1, int srid2, projPJ *pj1, projPJ *pj2) -{ - Proj4Cache *proj_cache = NULL; - - /* Set the search path if we haven't already */ - SetPROJ4LibPath(); - - /* get or initialize the cache for this round */ - proj_cache = GetPROJ4Cache(fcinfo); - if ( !proj_cache ) - return LW_FAILURE; - - /* Add the output srid to the cache if it's not already there */ - if (!IsInPROJ4Cache(proj_cache, srid1)) - AddToPROJ4Cache(proj_cache, srid1, srid2); - - /* Add the input srid to the cache if it's not already there */ - if (!IsInPROJ4Cache(proj_cache, srid2)) - AddToPROJ4Cache(proj_cache, srid2, srid1); - - /* Get the projections */ - *pj1 = GetProjectionFromPROJ4Cache(proj_cache, srid1); - *pj2 = GetProjectionFromPROJ4Cache(proj_cache, srid2); - - return LW_TRUE; -} /** diff --git a/regress/out_geography_expected b/regress/out_geography_expected index a469e7385..987ec57b9 100644 --- a/regress/out_geography_expected +++ b/regress/out_geography_expected @@ -16,7 +16,7 @@ gml_prefix_01|1,2 gml_prefix_02|1 2 gml_prefix_03|1,2 gml_prefix_04|1 2 -ERROR: Only SRID 4326 is currently supported in geography. +ERROR: GetProj4StringSPI: Cannot find SRID (10) in spatial_ref_sys kml_srid_02|0,1 kml_empty_geom| kml_precision_01|1,1 @@ -26,7 +26,7 @@ ERROR: Only KML 2 is supported ERROR: Only KML 2 is supported kml_prefix_01|1,2 kml_prefix_02|1,2 -ERROR: Only SRID 4326 is currently supported in geography. +ERROR: Only lon/lat coordinate systems are supported in geography. svg_empty_geom| svg_option_01|M 1 -1 L 4 -4 5 -7 svg_option_02|M 1 -1 l 3 -3 1 -3 @@ -44,8 +44,8 @@ geojson_crs_01|{"type":"Point","crs":{"type":"name","properties":{"name":"EPSG:4 geojson_crs_02|{"type":"Point","crs":{"type":"name","properties":{"name":"EPSG:4326"}},"coordinates":[1,1]} geojson_crs_03|{"type":"Point","crs":{"type":"name","properties":{"name":"urn:ogc:def:crs:EPSG::4326"}},"coordinates":[1,1]} geojson_crs_04|{"type":"Point","crs":{"type":"name","properties":{"name":"urn:ogc:def:crs:EPSG::4326"}},"coordinates":[1,1]} -ERROR: Only SRID 4326 is currently supported in geography. -ERROR: Only SRID 4326 is currently supported in geography. +ERROR: GetProj4StringSPI: Cannot find SRID (1) in spatial_ref_sys +ERROR: GetProj4StringSPI: Cannot find SRID (1) in spatial_ref_sys geojson_bbox_01|{"type":"LineString","coordinates":[[1,1],[2,2],[3,3],[4,4]]} geojson_bbox_02|{"type":"LineString","bbox":[1,1,4,4],"coordinates":[[1,1],[2,2],[3,3],[4,4]]} geojson_bbox_03|{"type":"LineString","crs":{"type":"name","properties":{"name":"EPSG:4326"}},"bbox":[1,1,4,4],"coordinates":[[1,1],[2,2],[3,3],[4,4]]}