[R-sig-Geo] convert point and map coordinates UTM to WSG

Thu Apr 2 16:06:54 CEST 2009

Dear Roger,

thank you for your help. There is something I don't quite understand.
If I compare distances between the converted coordinates (SP.utm.wgs)
with the coordinates (SP.wgs) that I use to compare with (in pairs),
the distances are more variant than I expected. These go from
approximately 15 meters to 192. I thought that the error would be more
or less constant. Is it normal?

> spDistsN1(coordinates(SP.wgs),coordinates(SP.utm.wgs[1,]),longlat=TRUE)[1]*1000
[1] 25.66325
> spDistsN1(coordinates(SP.wgs),coordinates(SP.utm.wgs[2,]),longlat=TRUE)[2]*1000
[1] 15.8736
> spDistsN1(coordinates(SP.wgs),coordinates(SP.utm.wgs[3,]),longlat=TRUE)[3]*1000
[1] 158.1777
> spDistsN1(coordinates(SP.wgs),coordinates(SP.utm.wgs[4,]),longlat=TRUE)[4]*1000
[1] 194.8797

Thank you,

Marc

2009/4/1 Roger Bivand <Roger.Bivand at nhh.no>
>
> On Wed, 1 Apr 2009, Marc Marí Dell'Olmo wrote:
>
>> Dear all,
>>
>> I have tried to transform coordinates from UTM Zone-31 ED-50 Spain Portugal
>> (EUR-D code in this web
>> http://earth-info.nga.mil/GandG/coordsys/onlinedatum/CountryEuropeTable.html#
>> EURD) to coordinates WSG84 ( google maps coordinates) with the
>> "ptransform" of the proj4 library and I can not do it right (first syntax).
>> I have tried also using the "spTransform" of the rgdal  library (second
>> syntax) and I have not been successful.
>>
>> I attach the syntax used, I got 4 points with coordinates ED ESP-50 ( "utm"
>> in the syntax) and their respective coordinates WSG84 (wgs.gold "in the
>> syntax). Can anyone tell me where is the error?
>
> So far several. One is reversing x and y values, another saying "latlong" when you mean "longlat", a third is the units of your input UTM (I think mm, not m), here, I'm going from longlat to UTM:
>
> library(rgdal)
>
> y.wgs <- c(41.4036923115558, 41.38010926480547, 41.41148169051035,
>  41.409534433320914)
> x.wgs <- c(2.145810127258301, 2.1807003021240234, 2.2183799743652344,
>  2.1170568466186523)
> SP <- SpatialPoints(cbind(x.wgs, y.wgs),
>  proj4string=CRS("+proj=longlat +datum=WGS84"))
> spTransform(SP,
>  CRS("+proj=utm  +ellps=intl +zone=31 +towgs84=-84,-107,-120,0,0,0,0,0"))
>
> The remaining differences are from the inadequacy of the +towgs84 string; Grids & Datums (http://www.asprs.org/resources/grids/) - like your source - gives +/- ranges for the transformation parameters, so you'll need to fish unless you can get a better authority (and there are few better than Grids & Datums).
>
> Hope this helps,
>
> Roger
>
>
>
>>
>> Thank you very much,
>>
>> Marc
>>
>>> #First sintax
>>> library(proj4)
>>> x.utm <- c(428727134,431587743,434639340,426268557)
>>> y.utm <- c(4584118048,4581465377,4584855001,4584607449)
>>> x.wgs <-
>>
>> c(41.4036923115558,41.38010926480547,41.41148169051035,41.409534433320914)
>>>
>>> y.wgs <-
>>
>> c(2.145810127258301,2.1807003021240234,2.2183799743652344,2.1170568466186523)
>>
>>>
>>> utm <- cbind(x.utm, y.utm)
>>> wgs.new <- ptransform(utm, "+proj=utm  +ellps=intl +zone=31
>>
>> +towgs84=-84,-107,-120,0,0,0,0,0 +no_defs",
>> +                  "+proj=latlong  +datum=WGS84")
>>>
>>> wgs.new
>>
>>         [,1]     [,2]     [,3]
>> [1,] -2.236353 1.570775 39.63333
>> [2,] -2.236353 1.570775 39.63333
>> [3,] -2.236353 1.570775 39.63333
>> [4,] -2.236353 1.570775 39.63333
>>>
>>> (wgs.gold <- cbind(x.wgs, y.wgs))
>>
>>       x.wgs    y.wgs
>> [1,] 41.40369 2.145810
>> [2,] 41.38011 2.180700
>> [3,] 41.41148 2.218380
>> [4,] 41.40953 2.117057
>>>
>>>
>>>
>>> #Second sintax
>>> library(rgdal)
>>> (data.utm <- data.frame(x.utm,y.utm,x.wgs,y.wgs))
>>
>>     x.utm      y.utm    x.wgs    y.wgs
>> 1 428727134 4584118048 41.40369 2.145810
>> 2 431587743 4581465377 41.38011 2.180700
>> 3 434639340 4584855001 41.41148 2.218380
>> 4 426268557 4584607449 41.40953 2.117057
>>>
>>> coordinates(data.utm) <- c("x.utm", "y.utm")
>>> proj4string(data.utm) <- CRS("+proj=utm  +ellps=intl +zone=31
>>
>> +towgs84=-84,-107,-120,0,0,0,0,0 +no_defs")
>>>
>>> (wgs.new2 <- spTransform(data.utm, CRS("+proj=latlong +datum=WGS84")))
>>
>>         coordinates    x.wgs    y.wgs
>> 1 (-128.134, 89.9988) 41.40369 2.145810
>> 2 (-128.134, 89.9988) 41.38011 2.180700
>> 3 (-128.134, 89.9988) 41.41148 2.218380
>> 4 (-128.134, 89.9988) 41.40953 2.117057
>>
>>
>>
>> 2009/3/23 Paul Hiemstra <p.hiemstra at geo.uu.nl>
>>
>>> Hi,
>>>
>>> To transform coordinates you can use the spTransform function from the
>>> rgdal pacakge. It uses proj4 to do the transformations. Proj uses a string
>>> to describe a certain projection, the so called proj4string. If you can
>>> determine the proj4string for both your projections you can transform
>>> between them. The first string (WGS84) looks like this:
>>>
>>> "+proj=longlat +datum=WGS84"
>>>
>>> The second (UTM) something like:
>>>
>>> "+proj=utm +zone=yourzone +datum=yourdatum"
>>>
>>> See the documentation of spTransform for more information on how to use the
>>> function.
>>>
>>> cheers and hth,
>>> Paul
>>>
>>>
>>> Marc Mar? Dell'Olmo wrote:
>>>
>>>> Dear R-users,
>>>>
>>>> Can anyone tell me how to convert point coordinates from WSG84 (google
>>>> maps coordinates) format to UTM HUSO-31 ED-50. And if I have a map
>>>> (shp format) with UTM HUSO-31 ED-50 coordinates, how can I convert to
>>>> WSG84 (google maps coordinates) coordinates. Should I install some
>>>> external software?
>>>>
>>>> Can anyone tell me the instructions?
>>>>
>>>> Thank you,
>>>>
>>>> Marc
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at stat.math.ethz.ch
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>
>>>>
>>>
>>>
>>> --
>>> Drs. Paul Hiemstra
>>> Department of Physical Geography
>>> Faculty of Geosciences
>>> University of Utrecht
>>> Heidelberglaan 2
>>> P.O. Box 80.115
>>> 3508 TC Utrecht
>>> Phone:  +3130 274 3113 Mon-Tue
>>> Phone:  +3130 253 5773 Wed-Fri
>>> http://intamap.geo.uu.nl/~paul <http://intamap.geo.uu.nl/%7Epaul>
>>>
>>>
>>
>>        [[alternative HTML version deleted]]
>>
>>
>
> --
> Roger Bivand
> Economic Geography Section, Department of Economics, Norwegian School of
> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: Roger.Bivand at nhh.no