Интерполировать нерегулярные точки данных x,y в регулярную сетку для отображения контуров
Я геолог, которому нужно создать несколько сотен согласованных контурных карт в проекте с различными наборами данных x y z.
Контурирование нерегулярных точек данных x y z включает создание "сетки" интерполированных (экстраполированных) значений z в однородной сетке xy. За пределами R - этот шаг называется "сеткой". Я относительно новичок в R и пытаюсь установить сильный рабочий процесс для построения сетки большого количества нерегулярных точек данных. Я борюсь!
В классическом программном обеспечении для контурных карт и в рабочем процессе выполняются следующие шаги:
- Читать данные XYZ
- Определите область интереса (AOI) для окончательной карты. XMIN, XMAX, YMIN, YMAX
- Определить сеточные интервалы (XINT, YINT) - устанавливает количество строк и количество столбцов в "сетку" (NROW, NCOL)
- Примените один из необходимых интерполяторов - которые создают 'z' на регулярной сетке / сетке (общие интерполяторы: обратное расстояние, обратное квадратное расстояние, средневзвешенное значение, полином, кригинг, сплайн и т. Д.)
- Контур получившейся "сетки"
Я пытаюсь написать скрипт R, чтобы точно следовать приведенной выше последовательности шагов для гибкости и контроля на протяжении всего анализа.
df - кадр данных, состоящий из примера набора данных.
wellid property z x y
060010 1 0.008849558 756994.5 2637732
009410 1 0.260162602 760190.9 2622262
009910 1 0.115044248 760898.7 2637466
051110 1 0.109243697 761690.2 2630985
065610 1 0.066666667 763064.1 2620929
011010 1 0.000000000 763089.3 2630888
035210 1 0.022556391 765942.4 2625944
052510 1 0.157894737 767058.1 2650034
006610 1 0.045045045 768265.0 2645318
009010 1 0.378151261 768471.8 2636731
011210 1 0.028776978 771393.8 2629001
064810 1 0.428571429 771394.1 2650776
009110 1 0.064220183 775332.6 2648531
011410 1 0.148760331 778324.8 2633905
065010 1 0.514851485 780480.9 2654874
052410 1 0.173913043 780961.0 2637571
064110 1 0.019417476 781001.5 2650994
009310 1 0.037383178 783904.7 2641130
010810 1 0.041237113 786200.6 2652417
052610 1 0.150537634 788007.5 2654005
Область интереса определяется из области исследования, как показано ниже:
xmin <- signif(min(wellcoords$x),4) - 1000
xmax <- signif(max(wellcoords$x),4) +1000
ymin <- signif(min(wellcoords$y),5) - 1000
ymax <- signif(max(wellcoords$y),5) +1000
xrange <- xmax-xmin
yrange <- ymax-ymin
gridint <- 500 # grid interval is set same for xint and yint
Значения: 754700, 791500,26196000,2658600, 36800, 39000, 500 соответственно.
После многих неудачных попыток - получил функцию interp() из пакета - akima для выполнения необходимой интерполяции. Благодаря ответу в разделе "Построение контуров на неправильной сетке"
fld<- with(df, interp(x=df$x, y=df$y, z=df$z, xo=xcoord, yo=ycoord, linear = FALSE, extrap = TRUE))
Это не позволило мне указать нужные элементы управления AOI. Я попытался использовать пакет MBA и все же работал над созданием параметра xy.est (сетка сетки) в качестве обязательного ввода.
Если генерируется правильная "сетка", ggplot2 и другие функции отображения являются мощными и достаточными.
Есть ли подходящие пакеты "Gridding" или "Шаги". Заранее спасибо.
1 ответ
Я не думаю, что вы должны использовать другой пакет, чем akima
(и графический пакет вроде ggplot2
). Вы можете указать AOI и No 'grid' в качестве аргументов interp. xo
а также yo
, И вы можете получить параметр xy.est interp2xyz(interp.obj)
,
df <- "your example data set"
# I didn't know What wellcoords were, so I treated df as wellcoords. These values are different from what you said.
xmin <- signif(min(df$x),4) - 1000 # 756000
xmax <- signif(max(df$x),4) + 1000 # 789000
ymin <- signif(min(df$y),5) - 1000 # 2619900
ymax <- signif(max(df$y),5) + 1000 # 2655900
gridint <- 500
library(akima)
fld<- with(df, interp(x = x, y = y, z = z, linear = FALSE, extrap = TRUE,
xo=seq(xmin, xmax, length=gridint),
yo=seq(ymin, ymax, length=gridint))) # give AOI and NO of 'grid'
# check whether the conditions are met.
length(fld$x); length(fld$y); length(fld$z); range(fld$x); range(fld$y)
# 500, 500, 250000 (=500^2), 756000 789000, 2619900 2655900, # all OK
contour(fld) # Left graph (most basic graphic output)
fld2 <- as.data.frame(interp2xyz(fld)) # the xy.est parameter (data.frame)
library(ggplot2)
ggplot(fld2, aes(x=x, y=y, z=z)) + geom_contour() # Right graph (simple example)