[R-sig-Geo] what Manifold does that Rgdal fails? Projection transformation R

[Roger Bivand]

On Mon, 22 May 2017, Slawomir Kozielec wrote:

> Dear friends,
>
> pretty new to the topic so let me know if i do not follow and reproducible
> examples convention. I am OK with R, but spatial data is still terra
> incognita

Please find out about ellipsoids and datums (sorry, not plural data). They 
describe the assumed shape of the world, and the relative distance of a 
point on the surface to an assumed 3D origin (roughly).

In one of your CRS, the ellipsoid is Airy, in another where a=b it is a 
sphere. In addition, EPSG 4269 presupposes NAD 1983, equivalent to GRS 
1980 and effectively WGS84. What you do not have are +towgs84= 
definitions, but if these are not given, they are assumed by proj in rgdal 
to be WGS84.

None of these are "correct", all are based on chosen assumptions. Why 
Manifold is assuming a spheroid is unknown. Define your datums correctly, 
and things should clear up.

Hope this clarifies,

Roger

>
> I got a Mapinfo file and i had to convert it to ESRI shp and then within r
> into spatial lines data frame and use with leaflet+OS maps epsg:3857
>
> whatever i tried within R with spTransform (rgdal) it did not work - either
> was not visible or was visible with offset (approx 2 meters). Proj4 did not
> want to work with SLDF.
> I took the file into Manifold, converted it to 3857, returned to R,
> transformed it to epsg:4269 and it works perfectly OK with leaflet + OS map
> 3857. But i have to automate this process to work within R, without any
> Manifold intervention.
>
> I checked how the files differed:
> i selected one item (the same of course) from the file and to make it easy
> to present created a centroid for it. Then I checked projection, xy, and
> latlong.
>
> the item created from the file without using Manifold:
> proj4string
>
>    "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000
> +y_0=-100000 +ellps=airy +units=m +no_defs"
>
> centroid:
>
>    528685 181836.2
>
> centroid after spTransform to 4269
>
>    -0.1450149 51.52028
>
> the item from the file pre-processed by Manifold:
> proj4string
>
>    "+proj=merc +lon_0=0 +lat_ts=0 +x_0=0 +y_0=0 +a=6378137 +b=6378137
> +units=m +no_defs"
>
> centroid:
>
>    -16321.64 6713938
>
> centroid after spTransform to 4269:
>
>    -0.1466198 51.52079
>
> the latter is correct of course. Just to emphasize - it is the same file
> and the same item
>
> Few questions:
>
> 1. any idea what happens in Manifold that fails in rgdal?
>
> 2. is there any mechanism in R that I can apply to fix the error (to match
> the 'with Manifold' output)? Maybe a kind of offset or work on prj files?
>
> I tried to use the latter CRS for transformation
>
>    P5 <-spTransform(ptp.points1P, CRS( "+proj=merc +lon_0=0 +lat_ts=0
> +x_0=0 +y_0=0 +a=6378137 +b=6378137 +units=m +no_defs"))
>
> and it returned still incorrect values
> the xy is t least inverted but shifted
>
>     -16142.98 6713847
>
> the latlong is of course incorrect, but the same as without any attempt to
> transform
>
>    -0.1450149 51.52028
>
>
> content of prj file not transformed - failing to correctly show
>
>    PROJCS["Transverse_Mercator",GEOGCS["GCS_Airy
> 1830",DATUM["D_unknown",SPHEROID["airy",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]]
>
> content of prj manifold transformed file - good to go
>
>    PROJCS["Mercator_2SP",GEOGCS["GCS_unnamed
> ellipse",DATUM["D_unknown",SPHEROID["Unknown",6378137,0]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Mercator"],PARAMETER["standard_parallel_1",0],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]
>
> thanks in advance
>
> 	[[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
>

-- 
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: Roger.Bivand at nhh.no
Editor-in-Chief of The R Journal, https://journal.r-project.org/index.html
http://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en