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

Michael Sumner mdsumner at gmail.com
Tue May 23 11:11:08 CEST 2017


Appreciate the feedback, cheers

(proj4 can be willed into use for decomposed sp/sf/etc but beware of
internal transpose assumption in current version, going rogue is generally
not advisable)


On Tue, 23 May 2017, 18:48 Slawomir Kozielec <kozielec at gmail.com> wrote:

> all, Thank you for your assistance,
> I resolved the issue, though I am not sure how to 'call off' my request
> for advice.
> Basically the MapInfo was projected to BNG (27700), so rgdal tried to
> transform already projected file. On Manifold I most likely and purely by
> chance set it to initial BNG and this is how it was produced correctly (i
> have no practical experience with Manifold - I got it once upon a time with
> intention to learn but never had time to). Probably Proj4 could have worked
> as well (as you provide the initial projection), but it refuses to work
> with spatial dataframes, especially lines (as points you can just extract
> as coords and transform).
> What I did was basically to replace by substitution (not transformation)
> the proj4file attribute to BNG and then transformed it - everything went
> smoothly from then on
>
> Mike, the file format is not my choice - I am preparing an application for
> my team and shp is our internal standard. Alas our GIS team leader took few
> months off, some project to clean oceans off plastic, and his deputees did
> not know how to arrange the file for me - so i took a 'raw' file from our
> vendor. The file was in TAB format. I am not disputing whatever is better
> as I am a proper layman in the order of datums, elipsums and all other
> Latin words I would love to know:)
>
>
> On Tue, May 23, 2017 at 8:40 AM, Michael Sumner <mdsumner at gmail.com>
> wrote:
>
>> Quickly, there is no need to convert to shapefile. MapInfo TAB or MIF are
>> both vastly superior to SHP and are well supported by rgdal and sf.
>>
>> That takes out one weak link. I'm happy to explore with Manifold if you
>> provide reproducible examples. Also I would ask on georeference.org too,
>> though you'll get a very disparaging (of open source) perspective from some
>> of their staff. It works both ways.
>>
>> Cheers, Mike
>>
>> On Tue, 23 May 2017, 17:13 Roger Bivand <Roger.Bivand at nhh.no> wrote:
>>
>>> 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 <+47%2055%2095%2093%2055>; 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
>>>
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>> --
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> Kingston Tasmania 7050 Australia
>>
>>
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list