Как буферизовать растровые пиксели по их значениям?

28

Пиксели слева представляют местоположения дерева и связанные с ними радиусы короны (то есть значения пикселей в диапазоне от 2 до 5). Я хотел бы буферизовать эти растровые пиксели по значению радиуса короны. Изображение справа - это то, чего я надеюсь достичь, используя только методы растровой обработки .

Сначала я бы подумал об использовании круговой фокусной суммы в ArcGIS, хотя настройка соседства является фиксированным значением, которое не учитывает радиус коронки переменного размера.

Что такое хороший метод для «буферизации» пикселей по их значениям?

введите описание изображения здесь

Аарон
источник
2
Вы пытались преобразовать растр в точки, затем буфер за полем, а затем преобразовать обратно в растр?
2
Это помогает понять, что это нелокальная операция, потому что эта нелокальность показывает, что существуют внутренние ограничения на то, как можно выполнить расчет. Например, ваш вывод радикально изменится почти везде, если только один изолированный пиксель на входе изменится на большое значение. Таким образом, если вам известны ограничения входных значений, поделитесь ими, поскольку это может привести к улучшению решений. Например, будут ли все ваши входные значения всегда в наборе {2,3,4}?
whuber
@Dan Patterson Вот как я придумал изображение справа. Однако я стараюсь вообще избегать векторных операций и избегать этих шагов.
Аарон
2
@whuber Этот набор данных представляет деревья с переменным диаметром кроны. Учитывая это, измерения радиуса кроны дерева могут реально варьироваться от 1 до 10. Я должен также добавить, что буферизованные выходные данные должны быть только 0 для отсутствия короны и 1 для присутствия короны.
Аарон
1
Хорошо спасибо. Похоже, вы создали пример вывода, объединив 3-буфера точек со значением 3, 4-буфера точек со значением 4 и 5-буферы точек со значением 5. (Вы, кажется, забыли для обработки точек со значением 2.) Этот процесс не только отвечает на ваш вопрос, но (я полагаю) это самое простое решение с использованием инструментов, доступных в Spatial Analyst.
whuber

Ответы:

14

Вот чистое растровое решение по Python 2.7использованию numpyи scipy:

import numpy as np
from scipy import ndimage
import matplotlib.pyplot as plt

#create tree location matrix with values indicating crown radius
A = np.zeros((120,320))
A[60,40] = 1
A[60,80] = 2
A[60,120] = 3
A[60,160] = 4
A[60,200] = 5
A[60,240] = 6
A[60,280] = 7

#plot tree locations
fig = plt.figure()
plt.imshow(A, interpolation='none')
plt.colorbar()

#find unique values
unique_vals = np.unique(A)
unique_vals = unique_vals[unique_vals > 0]

# create circular kernel
def createKernel(radius):
    kernel = np.zeros((2*radius+1, 2*radius+1))
    y,x = np.ogrid[-radius:radius+1, -radius:radius+1]
    mask = x**2 + y**2 <= radius**2
    kernel[mask] = 1
    return kernel

#apply binary dilation sequentially to each unique crown radius value 
C = np.zeros(A.shape).astype(bool)   
for k, radius in enumerate(unique_vals):  
    B = ndimage.morphology.binary_dilation(A == unique_vals[k], structure=createKernel(radius))
    C = C | B #combine masks

#plot resulting mask   
fig = plt.figure()
plt.imshow(C, interpolation='none')
plt.show()

Входные данные: введите описание изображения здесь

Выход: введите описание изображения здесь

jatobat
источник
1
+1 за подход к расширению! Это работает и с близкими точками тоже.
Антонио Фальчано
Это отличный пример того, почему этот старый реактивный цвет был ужасен. Это выглядит намного яснее с viridis.
naught101
8

Векторный подход

Эта задача может быть выполнена в три этапа:

Примечание: использование поля буфера позволяет избежать расчета буфера для каждого значения радиуса короны.


Растровый подход

Избегая векторного решения, эта проблема предлагает использовать своего рода сотовые автоматы на основе ближайших соседей. Предполагая, что все черные пиксели являются нулями, пиксели возводятся в квадрат и их размер равен 1 (или, альтернативно, масштабируется своевременно), правила для принятия очень просты:

  1. Если значение пикселя ( VALUE) больше 1, оно становится равным, VALUE-1а затем учитывает окружающие пиксели. Если их значения меньше чем VALUE-1, эти пиксели рождаются или растут, и их значение становится VALUE-1. В противном случае эти пиксели сохраняются и остаются без изменений.
  2. Если VALUE<=1ничего не делать (пиксель мертв!).

Эти правила должны применяться до тех пор, пока все пиксели не станут мертвыми, то есть их значения равны 0 или 1. То есть N-1, где Nмаксимальное значение, которое вы имеете во входном растре. Этот подход может быть довольно легко реализован с немного Python и Numpy.

Антонио Фальчано
источник
1
Спасибо за ответ afalciano. Этот метод - то, как я создал изображение справа и использует векторный подход - тот, которого я пытаюсь избежать.
Аарон
1
Хорошо, Аарон, вот подход, основанный на растре. Надеюсь это поможет.
Антонио Фальчано
7

Другой вариант - создать отдельные растры для каждого значения пикселя, в данном случае 4 растра, с условием. Затем разверните растры на количество пикселей, соответствующее значению растра (возможно, перебирая список значений). Наконец, объедините растры (алгебраические или пространственные), чтобы создать один двоичный растр для крон деревьев.

HDunn
источник
1
Эта идея правильная. Детали могут быть улучшены: (1) выбор создает двоичный (0,1) индикатор деревьев заданного радиуса кроны. (2) Фокальная сумма этого выбора - с использованием круговой окрестности заданного радиуса - быстро вычисляется с использованием БПФ. (3) Добавление фокусных сумм (по точкам) и сравнение их с 0 дает желаемый буфер.
whuber
7

Это сложный вопрос, чтобы сделать это в растре, потому что у вас нет возможности использовать значение пикселя для определения размера буфера. Поэтому вам нужно будет сделать фокусный фильтр для каждого значения, как вы уже сказали.

Вот возможный ответ, чтобы сделать это только с 3 фильтрами (я не мог найти меньше), но не совсем так, как упомянуто Whuber: ваши буферы будут усечены, когда деревья расположены близко друг к другу.

1) РЕДАКТИРОВАТЬ: Евклидово распределение (это не полностью решает проблему, так как оно сокращает буферы в окрестностях меньших деревьев, но это лучше, чем артефакты моего первого решения).

2) евклидово расстояние вокруг каждого пикселя

3) растровый калькулятор (алгебра карт) с условным оператором

Con("allocation_result" > ("distance_result" / pixel_size) , 1 , 0)

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

radouxju
источник
+1 Это творческий подход. Я проверю, чтобы увидеть, возможно ли увеличить масштаб с использованием этого подхода.
Аарон
2
Евклидово дистанционное приближение не будет работать, потому что оно вычисляет только расстояние до ближайшего дерева, которое не обязательно является расстоянием до дерева, крона которого покрывает точку.
whuber
2

Интересно, почему вы не используете инструмент расширения ArcGIS ?

import arcpy
from arcpy.sa import *

raster_in  = r'c:\test.tif'
raster_out = r'c:\test_out.tif'

outExpand1 = Expand(raster_in, 2, 2)
outExpand2 = Expand(outExpand1, 3, 3)
outExpand3 = Expand(outExpand3, 4, 4)
outExpand4 = Expand(outExpand4, 5, 5)

outExpand4.save(raster_out)

В случае перекрытия: последняя expandкоманда будет охватывать предыдущие.

Мистер Че
источник
2

Если у вас есть позиция в пикселях, радиус и алгоритм средней точки круга (вариант алгоритма Брезенхэма) дают подсказку. IMO, легко создать многоугольник из этого подхода, и я думаю, что это легко реализовать в Python. Объединение этого набора полигонов дает вам площадь покрытия.

huckfinn
источник
Я знаю, что вопрос не в этом, но вы хотите узнать больше о графических примитивах и заполнении полигонов линий сканирования? Для кругов это очень просто. Выпуклая комбинация - модное слово и так далее ...
huckfinn
Как это будет применяться с использованием базовых растровых операций?
whuber
Если вы попытаетесь обработать это в растровом пространстве, определите точки круга, отсортировав их по y или x, и заполните пространство прямой линией (линия сканирования), которая находится на пути к заполнению круга. В треугольном подходе, если вы строите круг с помощью аппроксимации триангулярных секторов и пытаетесь заполнить треугольник, вам нужен тест, если точка находится внутри или снаружи (выпуклая комбинация) и является другим способом. А в подходе «ГИС» создание многоугольников (ориентированных по часовой стрелке многоугольников) и создание союза является третьим (ИМО - самый дорогой компьютер).
huckfinn
Чтобы быть ясным: и в подходе «ГИС» ... сделать алгебраическую операцию, такую ​​как объединение, пересечение, прикосновение ... является третьим ИМО самым дорогим вычислительным методом.
huckfinn