извлекать значение растрового пикселя, а также координаты xy пикселя с помощью объекта SpatialLine в R

У меня есть Spatialline, который я преобразовал из шейп-файла многоугольника (оцифрованный вручную на основе функций в "imagebrick" - это означает, что пространственно "ломаная линия" и "imagebrick" перекрываются, как я хотел)

polyline <- as(shapefiles_data[1,],"SpatialLines")

> polyline
class       : SpatialLines 
features    : 1 
extent      : 357714.3, 357719, 4076030, 4076035  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=11 +datum=NAD27 +units=m +no_defs +ellps=clrk66 +nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat 

И растровый кирпич

> imagebrick
class       : RasterBrick 
dimensions  : 29180, 14940, 435949200, 4  (nrow, ncol, ncell, nlayers)
resolution  : 0.6, 0.6  (x, y)
extent      : 354038.4, 363002.4, 4066992, 4084500  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=11 +datum=NAD27 +units=m +no_defs +ellps=clrk66 +nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat 
names       : band1, band2, band3, band4 
min values  :     0,     0,     0,     0 
max values  : 65535, 65535, 65535, 65535 

Я пытался извлечь пиксели в Rasterbrick, которых коснулась линия (SpatialLine). Но я также хотел извлечь координаты xy, связанные с этими извлеченными пикселями, и сохранить их вместе в SpatialPointsDataFrame. Я пробовал:

extract_polyline_test <- extract(imagebrick, polyline)
extract_polyline_test <- as.data.frame(extract_polyline_test)
> extract_polyline_test
band1 band2 band3 band4
1    125   156    73   392
2    129   161    80   366
3    156   208   122   374
4    142   175    97   330
5    102   117    31   389
6    117   143    57   381
7    162   221   127   413
8    213   302   204   405
..   ...   ...   ...   ...

Затем я пытался использовать функцию RasterToPoints, чтобы вытащить координаты xy. Для этого требуется растровый слой, содержащий все пиксели, из которых я хочу xy. Я пытался использовать «кадрирование», чтобы получить эту растровую «линию», но это не сработало.

rastercrop_polyline <-crop(imagebrick,polyline)
rastercrop_polyline

> rastercrop_polyline
class       : RasterBrick 
dimensions  : 8, 7, 56, 4  (nrow, ncol, ncell, nlayers)
resolution  : 0.6, 0.6  (x, y)
extent      : 357714.6, 357718.8, 4076030, 4076035  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=11 +datum=NAD27 +units=m +no_defs +ellps=clrk66 +nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat 
data source : in memory
names       : band1, band2, band3, band4 
min values  :    72,    63,     0,   156 
max values  :   274,   415,   299,   497 

s.polyline <- rasterToPoints(rastercrop_polyline)
> s.polyline
         x       y band1 band2 band3 band4
 [1,] 357714.9 4076035    93    96    17   342
 [2,] 357715.5 4076035    72    63     0   377
 [3,] 357716.1 4076035    90    93    17   371
 [4,] 357716.7 4076035   125   156    73   392
 [5,] 357717.3 4076035   129   161    80   366
 [6,] 357717.9 4076035   156   208   122   374
 [7,] 357718.5 4076035   142   175    97   330
 [8,] 357714.9 4076034   100   108    25   326
 ..  ........  .......   ...   ..     ..    ...

Кажется, что функция «кадрирования» обрезает не «линию» растра, а целый «многоугольник». Есть ли другой способ получить координаты xy пикселей, которые я извлек?

p.s. Я также пытался извлечь окружающие пиксели из каждого извлеченного пикселя. Как мне это сделать? Я заметил, что в функции извлечения есть параметр функции. Должен ли я работать над написанием небольшой функции внутри функции извлечения или есть другие способы сделать это?

p.p.s Я также думал использовать функцию Rasterize для растеризации «SpatialLine», а затем использовать RasterToPoints для извлечения координат xy. Но я не могу контролировать, чтобы каждый пиксель «растеризованной SpatialLine» был равен пикселям, которых коснулась «линия» в исходном блоке изображения.

Спасибо всем за то, что уделили время просмотру моего вопроса, и заранее спасибо за вашу помощь.


person nomoreraster    schedule 18.03.2014    source источник


Ответы (1)


Вы говорите, что «думали использовать функцию Rasterize для растрирования SpatialLine, а затем использовать RasterToPoints для получения координат xy».

Это кажется правильным для того, что вы хотите

Но затем вы говорите, что «не можете контролировать каждый пиксель« растеризованной SpatialLine », равный пикселям, которых касается« линия »в исходном блоке изображения».

Но они есть, за исключением того, что ячейки с несколькими линиями представлены только один раз.

Другой подход:

# set up some example data
r <- raster()
values(r) <- runif(ncell(r))
cds1 <- rbind(c(-50,0), c(0,60), c(40,5), c(15,-45), c(-10,-25))
cds2 <- rbind(c(80,20), c(140,60), c(160,0), c(140,-55))
lines <- SpatialLines(list(Lines(list(Line(cds1)), "1"), Lines(list(Line(cds2)), "2") ))

# values for lines
v <- extract(r, lines)
# create raster with cell numbers
ids <- init(r, v='cell')
# extract these
cells <- extract(ids, lines)
# compute xy
xy = lapply(cells, function(x) xyFromCell(r, x))
person Robert Hijmans    schedule 18.03.2014
comment
Роберт, спасибо за помощь! Я реализовал ваш подход, и он хорошо сработал. Большое тебе спасибо! - person nomoreraster; 19.03.2014
comment
Кроме того (как Роберт, конечно, уже знает;), можно получить объединенные значения v и cells за один шаг с помощью vv <- extract(r, lines, cellnumbers=TRUE). - person Josh O'Brien; 06.03.2015