what Manifold does that Rgdal fails? Projection transformation R

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
5 messages Options
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

what Manifold does that Rgdal fails? Projection transformation R

Slawomir Kozielec
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

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
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: what Manifold does that Rgdal fails? Projection transformation R

Roger Bivand
Administrator
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
> [hidden email]
> 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: [hidden email]
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
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand
Department of Economics
Norwegian School of Economics
Helleveien 30
N-5045 Bergen, Norway
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: what Manifold does that Rgdal fails? Projection transformation R

Michael Sumner-2
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 <[hidden email]> 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
> > [hidden email]
> > 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: [hidden email]
> 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
> [hidden email]
> 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

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: what Manifold does that Rgdal fails? Projection transformation R

Slawomir Kozielec
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 <[hidden email]> 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 <[hidden email]> 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
>> > [hidden email]
>> > 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:
>> [hidden email]
>> 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
>> [hidden email]
>> 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
>
>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Threaded
Open this post in threaded view
|  
Report Content as Inappropriate

Re: what Manifold does that Rgdal fails? Projection transformation R

Michael Sumner-2
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 <[hidden email]> 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 <[hidden email]>
> 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 <[hidden email]> 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
>>> > [hidden email]
>>> > 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:
>>> [hidden email]
>>> 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
>>> [hidden email]
>>> 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]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Loading...