[postgis-users] Implementation of postgis interoperability with R

Paul Ramsey pramsey at cleverelephant.ca
Thu Mar 12 19:36:23 PDT 2020


proj_identify:

https://proj.org/development/reference/functions.html#_CPPv413proj_identifyP10PJ_CONTEXTPK2PJPKcPPCKcPPi

On Thu, Mar 12, 2020 at 7:35 PM Paul Ramsey <pramsey at cleverelephant.ca> wrote:
>
> (a) do you already know the authority and code? Then first use those
> to look-up the SRID in spatial_ref_sys (select srid from
> spatial_ref_sys where auth_name = 'EPSG' and auth_code = 14234)
> (b) do you only have WKT input? Then, sInce you have access to proj on
> the R side, I'd suggest instantiating your const PJ *obj, and then
> using proj_identify() to see if proj can look up your projection with
> a good match. Maybe it can!
> (c) try a text lookup, it will probably fail, since as you note
> there's a million ways text representations can differ while still
> being the same.
>
> P.
>
> On Thu, Mar 12, 2020 at 8:20 AM Etienne B. Racine <etiennebr at gmail.com> wrote:
> >
> > Dear list,
> >
> > I'm currently updating the simple feature package for R to use the latest proj6 to interact with postgis. One of the things we perform behind the scene is synching the projections between R and postgis.
> >
> > The current version of the package provides interoperability in two ways: read and write data to postgis, and allow the user to write R code that is executed, as SQL, on the remote database, just like a normal data.frame (R tables).
> >
> > Our current implementation relies on epsg code and proj4string matching to synch database srid with the local crs, but with PROJ6, we need to add proj wkt matching.
> >
> > Two questions:
> >
> > This is how I plan to adapt the projection synching, any suggestion for improvement or pointer to other implementations would be much appreciated:
> > 1. Match local epsg code and srid, create a crs in R using the `srtext` and check that the two projections are equivalent using `IsSame` from OGRSpatialReference.
> > 2. No matching code, so check `spatial_ref_sys` table to see if a proj WKT would match any `srtext` (slow and not very robust, I have a follow up question)
> > 3. No matching code and no WKT match, so try to match proj4string.
> > 4. Clearly no matching projection, so insert new srid in `spatial_ref_sys` and set `srid = max(srid) + 1`, `auth_name = 'sf'` (the simple feature package in R is called sf).
> > 5. If it fails (e.g. lack of permission), then error and ask the user to change the projection, or use srid = 0.
> >
> > Is string matching `=` the best way to match a `srtext` and `proj4string` in  `spatial_ref_sys`?
> > ``` sql
> > select * from spatial_ref_sys where srtext = {wkt};
> > ```
> > But `=` obviously rejects matches if parameters are ordered differently (and I suspect we could have issues with encoding). One solution is to read `spatial_ref_sys` to R, and match there, but it is very very slow since it requires to parse all the wkt to create a projection, and then perform a comparison with OGRSpatialReference->IsSame on every projection.
> >
> > Thanks for your help!
> > Etienne
> > _______________________________________________
> > postgis-users mailing list
> > postgis-users at lists.osgeo.org
> > https://lists.osgeo.org/mailman/listinfo/postgis-users


More information about the postgis-users mailing list