2015-05-24 6 views
0

Мои данные представлены в виде широты и долготы в виде градусов, минут и секунд. Это было сохранено в текстовой форме (character-class, Игнорировать row.names).Как преобразовать долготу/широту в UTM? R. rgdal, sp

> pasaporte 
     Latitud  Longitud 
4 13°50¨52" sur 73°45¨12" oeste 
36 13°01¨ sur 75°05¨ oeste 
46 13°09¨26" sur 74°13¨22" oeste 

С помощью answer, представленной в списке рассылки, я преобразовал данные в десятичную форму ...

> pasaporte 
    Latitud Longitud 
4 13.84778 73.75333 
36 13.01667 75.08333 
46 13.15722 74.22278 

После этого я преобразовал data.frame к а SpatialPoints объекта ,

xy <- data.frame(cbind(round(pasaporte[, 'Longitud'], 5), round(pasaporte[, 'Latitud'], 5))) #Rounded the decimals out of doubt of it interfering later 
xy <- SpatialPoints(coords = xy, 
       proj4string = CRS('+proj=longlat +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs')) 

Затем приступил к преобразованию проекции через spTransform

xy_utm <- spTransform(xy, CRS('+proj=utm +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs +towgs84:0,0,0')) 

Я преобразовал первую координату в примере, используя this website (делает Shure, чтобы обеспечить отрицательный, так как координаты относятся к Югу и Западу) и предоставляет следующие координаты: 633726.46 mE 8468757.51 mN Зона 18L, что является правильным. В отличии от этого sp объекта после преобразования:

SpatialPoints: 
      X1  X2 
[1,] 12051333 14804033 
[2,] 12569894 14792671 
[3,] 12301330 14684870 
Coordinate Reference System (CRS) arguments: +proj=utm +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs +towgs84=0,0,0 

spTransform руководства по конкретно упоминается, что метаданные должны быть правильно предусмотрены. Мои данные используют эллипсоид WGS84 и находятся в зоне 18L. Точнее, epsg:32718. Где я ошибаюсь?

+0

Это следует задать здесь: http://gis.stackexchange.com/questions/tagged/r –

ответ

2

Я бы рекомендовал «проецировать» (то есть умножить на -1) точки в границах -180.0000, -90.0000, 180.0000, 90.0000, поскольку они лежат в южном и западном полушариях.

library(sp) 

pasaporte <- structure(list(Latitud = c(13.84778, 13.01667, 13.15722), 
          Longitud = c(73.75333, 75.08333, 74.22278)), 
         .Names = c("Latitud", "Longitud"), class = "data.frame", 
         row.names = c("4", "36", "46")) 

pasaporte <- pasaporte * -1 

pasaporte 
#  Latitud Longitud 
# 4 -13.84778 -73.75333 
# 36 -13.01667 -75.08333 
# 46 -13.15722 -74.22278 

# Conversion into SpatialPoints 
coordinates(pasaporte) <- ~Longitud+Latitud 
# Setting default projection 
proj4string(pasaporte) <- CRS('+init=epsg:4326') 

pasaporte 
# SpatialPoints: 
#  Longitud Latitud 
# 4 -73.75333 -13.84778 
# 36 -75.08333 -13.01667 
# 46 -74.22278 -13.15722 
# Coordinate Reference System (CRS) arguments: +init=epsg:4326 
# +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 

# Projection 
# Be sure to have rgdal first installed. 
xy_utm <- spTransform(pasaporte, CRS('+init=epsg:32718')) 
xy_utm 
# SpatialPoints: 
# Longitud Latitud 
# 4 634726.5 8468758 
# 36 490964.2 8561019 
# 46 584231.8 8545348 
# Coordinate Reference System (CRS) arguments: +init=epsg:32718 +proj=utm 
# +zone=18 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 
# +towgs84=0,0,0 
+0

correct; это будет работать, однако, если установлен пакет 'rgdal'. –

+0

@ EdzerPebesma Точно. Спасибо за указание на это. –

Смежные вопросы