[R-sig-Geo] WKT coordinate system string (PRJ) to PROJ.4

[Roger Bivand]

On Tue, 5 Jan 2016, boB Rudis wrote:

> rgdal 1-1.4 (built from rforge) has "showP4" (and "showWKT"):
>
>  prj_1 <- 'PROJCS["Transverse_Mercator",
>    GEOGCS["GCS_OSGB 1936",
>    DATUM["D_OSGB_1936",
>    SPHEROID["Airy_1830",6377563.396,299.3249646]],
>    PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],
>    PROJECTION["Transverse_Mercator"],
>    PARAMETER["latitude_of_origin",49],
>    PARAMETER["central_meridian",-2],
>    PARAMETER["scale_factor",0.9996012717],
>    PARAMETER["false_easting",400000],
>    PARAMETER["false_northing",-100000],
>    UNIT["Meter",1]]'
>
>  showP4(prj_1)
>  ## [1] "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000
> +y_0=-100000 +datum=OSGB36 +units=m +no_defs "

Thanks for checking. The functions do use morphFromESRI() if requested - 
I'll revisit this so that morphing from on import and to on export are 
settable options for relevant functions, even when there are no arguments. 
I find

http://www.gdal.org/classOGRSpatialReference.html#af91af2639702e2793daf93ffe231b577

http://www.gdal.org/classOGRSpatialReference.html#ad556dfdc04d9ec5f1714fc6b5e0eb6a6

etc. confusing, given the need to set GDAL_FIX_ESRI_WKT in some cases. 
Adding towgs84 values does need attention - contributions welcome!

Roger

>
> On Tue, Jan 5, 2016 at 6:20 AM, boB Rudis <bob at rudis.net> wrote:
>> It's too early I guess as it didn't occur to me that if writing it out
>> locally was naff using an external service is more so. The
>> SpatialReference "importFromESRI()" function doesn't seem to be mapped
>> in rgdal. That would be the function needed to avoid a
>> write-out/external call.
>>
>> On Tue, Jan 5, 2016 at 5:57 AM, boB Rudis <bob at rudis.net> wrote:
>>> This may not be optimal as it uses an external service:
>>>
>>>   prj_to_epsg <- function(prj) {
>>>
>>>     require(sp)
>>>     require(httr)
>>>     require(jsonlite)
>>>
>>>     res <- GET("http://prj2epsg.org/search.json",
>>>                query=list(exact=TRUE,
>>>                           error=TRUE,
>>>                           mode="wkt",
>>>                           terms=prj))
>>>
>>>     # one shld prbly do more error checking than this
>>>     stop_for_status(res)
>>>
>>>     dat <- fromJSON(content(res, as="text", flatten=TRUE))
>>>
>>>     # NOTE: there could be more in dat$codes if prj was ambiguous
>>>     CRS(paste0("+init=epsg:", dat$codes[1, "code"]))
>>>
>>>   }
>>>
>>>   prj <- paste0(readLines("110m_admin_1_states_provinces_shp.prj"))
>>>   prj_to_epsg(prj)
>>>   ## CRS arguments:
>>>   ##  +init=epsg:4326 +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84
>>>   ## +towgs84=0,0,0
>>>
>>>   prj_1 <- 'PROJCS["Transverse_Mercator",
>>>     GEOGCS["GCS_OSGB 1936",
>>>     DATUM["D_OSGB_1936",
>>>     SPHEROID["Airy_1830",6377563.396,299.3249646]],
>>>     PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],
>>>     PROJECTION["Transverse_Mercator"],
>>>     PARAMETER["latitude_of_origin",49],
>>>     PARAMETER["central_meridian",-2],
>>>     PARAMETER["scale_factor",0.9996012717],
>>>     PARAMETER["false_easting",400000],
>>>     PARAMETER["false_northing",-100000],
>>>     UNIT["Meter",1]]'
>>>   prj_to_epsg(prj_1)
>>>   ## CRS arguments:
>>>   ##  +init=epsg:27700 +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717
>>>   ## +x_0=400000 +y_0=-100000 +datum=OSGB36 +units=m +no_defs +ellps=airy
>>>   ## +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894
>>>
>>> On Tue, Jan 5, 2016 at 5:09 AM, Michael Sumner <mdsumner at gmail.com> wrote:
>>>> Hello, is there any R-level way to convert a WKT/PRJ *coordinate system*
>>>> string to PROJ.4 used by sp/rgdal?
>>>>
>>>> (Previously I've written out to SHP/PRJ, and read back in with rgdal, but
>>>> that's obviously a bit naff.)
>>>>
>>>> (I know you can rgeos::readWKT and wkb::readWKB for *geometry*).
>>>>
>>>> Cheers, Mike.
>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> Dr. Michael Sumner
>>>> Software and Database Engineer
>>>> Australian Antarctic Division
>>>> 203 Channel Highway
>>>> Kingston Tasmania 7050 Australia
>>>>
>>>>         [[alternative HTML version deleted]]
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at r-project.org
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 91 00
e-mail: Roger.Bivand at nhh.no