Конвертировать координаты NZMG в широту / долготу

Question

У меня есть набор координат NZ Map Grid, которые я хочу преобразовать в широту / долготу. Исходя из этого вопроса , вот что я пытался.

library(sp)
options(digits = 11) # to display to greater d.p.

Попытка 1:

proj4string <- "+proj=nzmg +lat_0=-41.0 +lon_0=173.0 +x_0=2510000.0 
+y_0=6023150.0 +ellps=intl +units=m"
p <- proj4::project(c(2373200, 5718800), proj = proj4string, inverse=T)

Попытка 2

dat <- data.frame(id = c(1), x = c(2373200) , y = c(5718800))
sp::coordinates(dat) = ~x+y

sp::proj4string(dat) = CRS('+init=epsg:27200') 
data_wgs84 <- spTransform(dat, CRS('+init=epsg:4326'))
print(data_wgs84)

ЕслиЯ запускаю свои координаты через инструмент преобразования координат Линца Я получаю немного другой результат, который является "истинным" результатом.

Results:
171.30179199  -43.72743909  # attempt 1 - ~200m off linz 
171.30190004, -43.72577765  # attempt 2 - a few meters off linz
171.30189464, -43.72576664  # linz

Основано на ответе Майка Т . Я должен использовать «метод преобразования сетки искажений», а он ссылается на «файл сдвига сетки nzgd2kgrid0005.gsb».

Мой вопрос: возможно ли сделать это преобразование, используя R без загрузки дополнительных файлов (nzgd2kgrid0005.gsb)? Я хочу поделиться своим кодом с другими без необходимости загружать дополнительные файлы.

Любой совет высоко ценится.

Answer 1

Оказывается, это довольно просто, если у вас установлен пакет rgdal, необходимый файл nzgd2kgrid0005.gsb включен, и вам не нужно ничего загружать.

Вам просто нужно использовать полную строку PROJ.4, как указано в ответе Mike T .

dat <- data.frame(id = c(1), x = c(2373200) , y = c(5718800))
sp::coordinates(dat) = ~x+y

proj4string <- "+proj=nzmg +lat_0=-41 +lon_0=173 +x_0=2510000 +y_0=6023150 
+ellps=intl +datum=nzgd49 +units=m +towgs84=59.47,-5.04,187.44,0.47,-0.1,1.024,-4.5993 
+nadgrids=nzgd2kgrid0005.gsb +no_defs"

sp::proj4string(dat) = sp::CRS(proj4string) 
data_wgs84 <- sp::spTransform(dat, sp::CRS('+init=epsg:4326'))
as.data.frame(data_wgs84)

id          x            y
1 171.3018946 -43.72576664

То же, что и результатLINZ инструмент для преобразования координат. Надеюсь, это сэкономит кому-то еще немного времени.