Проекция на несколько зон UTM в R - PullRequest
Новая
       34

Проекция на несколько зон UTM в R

0 голосов
/ 18 июня 2020

У меня есть набор точек GPS в десятичных градусах для 120 образцов из 13 точек отбора проб (данные ниже - центроиды) 

data.frame(x=c(122.4978, 122.1819, 113.7757, 123.0335, 113.0614, 116.8117, 113.6997, 118.8820, 116.8047, 114.0167, 115.0560, 118.6114, 113.3973), y=c(-16.90613, -17.89392, -23.10492, -16.57618, -25.63859, -20.56915, -25.76536, -19.02967, -19.79827, -21.91585, -21.62384, -20.28056, -26.03371))

Теперь мне нужно вычислить кратчайшее расстояние между каждой парой точек ( в км). Эта часть ясна, и я успешно перепроецировал точки в другом проекте, где все точки находились в одной зоне UTM, а затем использовал gdistance::shortestPath с растром и матрицей перехода. Эти точки здесь, однако, охватывают 3 зоны UTM (49 - 51).

Я, к сожалению, не очень опытен в этих вопросах, поэтому был бы признателен, если бы кто-нибудь указал мне на ресурс, который, в идеале, объясняет как сделать проекцию на несколько зон шаг за шагом. Или, если это неправильный подход, дайте мне подсказку, как действовать дальше.

1 Ответ

0 голосов
/ 18 июня 2020

Из вашего вопроса я не уверен, что ясно, а что нет. Но я просто покажу вам, как получить «кратчайшее расстояние между каждой парой точек (в км)»

Данные вашего примера 

df <- data.frame(
  x=c(122.4978, 122.1819, 113.7757, 123.0335, 113.0614, 116.8117, 113.6997, 118.8820, 116.8047, 114.0167, 115.0560, 118.6114, 113.3973),
  y=c(-16.90613, -17.89392, -23.10492, -16.57618, -25.63859, -20.56915, -25.76536, -19.02967, -19.79827, -21.91585, -21.62384, -20.28056, -26.03371)
)

Решение 

library(raster)
p <- pointDistance(df, lonlat=T)
p[1:3, 1:3] / 1000
#          [,1]     [,2] [,3]
#[1,]    0.0000       NA   NA
#[2,]  114.3593    0.000   NA
#[3,] 1141.3329 1049.202    0

Или эквивалент, но с заполненными обоими треугольниками матрицы. 

library(geosphere)
d <- distm(df)
d[1:3, 1:3] / 1000
#          [,1]      [,2]     [,3]
#[1,]    0.0000  114.3593 1141.333
#[2,]  114.3593    0.0000 1049.202
#[3,] 1141.3329 1049.2017    0.000

Еще один, с terra, возвращающий матрицу расстояний (dist объект). 

library(terra)
v <- vect(df, crs="+proj=longlat +datum=WGS84")
s <- distance(v)
class(s)
#[1] "dist"
s[1:3]/1000
#[1]  114.35935 1141.33295   67.79552
...