Извлечение R-растра по SpatialLine: соотнесите извлеченные значения с фактическим расстоянием

При извлечении значений растра вдоль SpatialLine в R, как связать эти значения с фактическим расстоянием вдоль этой линии?

Предположим, я хочу извлечь значение логотипа R в следующей строке:

library(raster)
r <- raster(system.file("external/rlogo.grd", package="raster"))
x=c(5, 95)
y=c(20, 50)
line = SpatialLines(list(Lines(Line(cbind(x,y)), ID="a")))
plot(r)
plot(line, add=TRUE)

Rlogo

Я могу извлечь значения и построить их график, но как заменить значения x (1:length(vals) ниже) на фактическое расстояние (начиная, например, с 0 с левой стороны строки)?

vals <- extract(r, line)[[1]]
plot(1:length(vals), vals, type='o')

Rlogo_val

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


r r-raster extraction
person ztl    schedule 06.03.2015    source источник

Ответы (2)


arrow_upward
0
arrow_downward

Я не уверен, что именно вы спрашиваете, но если вы ищете расстояния между самой левой координатой отрезка линии и центрами ячеек, через которые проходит линия, то вы можете найти расстояния следующим образом:

x <- extract(r, l, cellnumbers=TRUE)[[1]]
xy <- xyFromCell(r, x[,1]) # get cell coordinates where the line passes
start <- xy[which.min(xy[,1]),] # leftmost coordinate of the line
d <- apply(xy, 1, function(x, start) sqrt(sum((x-start)^2)), start=start) # find distances between the line segment start and the cells
plot(1:length(d), d, type='o')
person jvj    schedule 07.03.2015
comment
Спасибо @jvj. Предполагая, что извлеченные значения связаны с самой строкой, я больше рассматривал расстояния на самой строке (я знаю, что это не совсем правильно, а также зависит от method, используемого для извлечения, но не должно быть слишком большой проблемой, если разрешение сетки высокое по сравнению с длиной линии). Интуитивно я больше думал о расстоянии между самой левой координатой линии и проекциями центров ячеек на линии (я полагаю, это легко существует где-то еще в R) - но, возможно, в более эффективный / строгий / встроенный способ ... - person ztl; 09.03.2015
comment
Кстати, @jvj, разве ваш start выше НЕ (или, скорее всего, нет, за исключением случая, когда линия проходит через центр ячейки) на линии? Кроме того, не гарантируется, что start определяется как самая левая и самая нижняя ячейка, не так ли? Или я где-то не прав? - person ztl; 10.03.2015
comment
Чтобы получить расстояния на самой линии, вы, возможно, можете использовать geosphere :: dist2Line - person Robert Hijmans; 14.03.2015

arrow_upward
0
arrow_downward

Вот решение (частично на основе ввода @ jvj) путем попытки вычислить ортогональные проекции центров ячеек, предоставленных raster::extract на линии, а затем вычислить расстояния вдоль линии.

(Это сценарий для начинающих R, который, вероятно, легко улучшить, но, похоже, работает (и, конечно, только для растров с проекцией с учетом расстояний))

vals <- extract(r, line, cellnumbers=TRUE)[[1]]
cellsxy <- xyFromCell(r, vals[,1]) # coordinates of intersected cells (likely not ON the line)
linexy = spsample(line, 1000, "regular") # get the line as points
linexy <- matrix(cbind(linexy$x, linexy$y), ncol=2) # easier than Spatial object for later
orthoproj <- c() # to store the orthogonal projections of cells centres on the line
for (i in 1:nrow(cellsxy)) {
  xypt = cellsxy[i,]
  min.index <- which.min(spDistsN1(linexy, xypt)) 
  orthopt <- linexy[min.index, ] # orthogonal projections = smaller distance 
  orthoproj <- c(orthoproj, c(orthopt[1], orthopt[2]))
}
orthoproj <- matrix(orthoproj, ncol=2, byrow=T)
orthoproj <- data.frame(x=orthoproj[,1], y=orthoproj[,2])
orthoproj <- orthoproj[order(orthoproj[,1]),] # reorder with increasing distance
orthoproj <- data.frame(x=orthoproj$x, y=orthoproj$y)
start <- linexy[which.min(linexy[,1]),] # leftmost coordinate of the line
dists <- apply(orthoproj, 1, 
               function(xy, start) sqrt(sum((xy-start)^2)), 
               start=start) # distances between 'start' and  the orthogonal projections

plot(dists, rev(vals[,2]), type='o') # !! beware: order of 'vals' and 'dists' 
                    # depending on the order in which cellnumbers are returned 
                    # in raster::extract and the shape of your line !!

Rlogo_extractedValuesAlongLineWithDistance

person ztl    schedule 13.03.2015