Создать массив numPy из растрового изображения
Я пытаюсь превратить 4-полосное (RGB и инфракрасное) растровое изображение в массив numPy в ArcMap. После успешного преобразования в массив Numpy я хочу подсчитать количество пикселей, которые не имеют данных на изображении. При проверке в ArcMap эти пиксельные цвета помечаются как "Нет", они выглядят черными, но в них отсутствуют данные канала красного, зеленого и / или синего из диапазонов 1,2 или 3. Мне нужно их найти.
Вот что у меня так далеко:
import numpy
import os
myDir = "C:\\Temp\\temp"
# myFile = "4_pixel_test.tif"
myFile = "4band.tif"
# import 4band (R,G,B & nr Infrared) image
fName = os.path.join(myDir, myFile)
head, tail = os.path.split(fName)
# Convert Raster to Array, Default using LowerLeft as Origin
rasArray = arcpy.RasterToNumPyArray(fName)
# find out the number of bands in the image
nbands = rasArray.shape[0] # int
# print nbands (int)
blackCount = 0 # count the black pixels per image
th = 0 # Threhold value
# print rasArray
r, g, b, a = rasArray # not working
rCheck = numpy.any(r <= th)
gCheck = numpy.any(g <= th)
bCheck = numpy.any(b <= th)
aCheck = numpy.any(a == 0)
print rCheck
print gCheck
print bCheck
print aCheck
# show the results
if rCheck:
print ("Black pixel (red): %s" % (tail))
elif gCheck:
print ("Black pixel (green): %s" % (tail))
elif bCheck:
print ("Black pixel (blue): %s" % (tail))
else:
print ("%s okay" % (tail))
if aCheck:
print ("Transparent pixel: %s" % (tail))
Ошибка во время выполнения (последний вызов был последним): файл "", строка 14, в файле "c:\program files (x86)\arcgis\desktop10.2\arcpy\arcpy__init__. Py", строка 1814, в RasterToNumPyArray возвращает _RasterToNumPyArray(*args, **kwargs) RuntimeError: ОШИБКА 999998: непредвиденная ошибка.
# previous code which might have incorrect numpy import
# options so I'm going with default options until I know better
# import numpy
# import os
#
# myDir = "C:\\Temp\\temp"
# myFile = "4_pixel_test.tif"
# fName = os.path.join(myDir, myFile)
#
# Convert Raster to Array
# rasArray = arcpy.RasterToNumPyArray(fName)
# maxVal = rasArray.max()
# minVal = rasArray.min()
# maxValpos = numpy.unravel_index(rasArray.argmax(),rasArray.shape)
# minValpos = numpy.unravel_index(rasArray.argmin(),rasArray.shape)
#
# desc = arcpy.Describe(fName)
# utmX = desc.extent.upperLeft.X + maxValpos[0]
# utmY = desc.extent.upperLeft.Y - maxValpos[1]
#
# for pixel in numpy.nditer(rasArray):
# # r,g,b = pixel # doesn't work - single dimension array
# print pixel
#
Я смог изменить растровое изображение в массив numPY из кода здесь.
Не уверен, как хранится массив numPY, но при итерации по нему данные распечатываются, начиная с оси y и работая вниз (столбец за столбцом) изображения вместо x (строка за строкой).
Мне нужно переключить это, чтобы я мог читать данные попиксельно (RGBA) сверху вниз, слева направо. Тем не менее, я не знаю достаточно о NumPy, чтобы сделать это.
Я думаю, что рассматриваемая ошибка может быть из-за размера рассматриваемого tiff: он отлично работает с 2,5 МБ, но падает на 4 ГБ.:(
2 ответа
Кажется, ты спрашиваешь о np.nditer,
Вы не хотите использовать nditer если вам не нужен контроль низкого уровня. Тем не менее, вам почти никогда не понадобится такой уровень контроля. Лучше не использовать nditer если вы не знаете точно, зачем вам это нужно.
То, что у вас есть, это трехмерный массив. В настоящее время вы перебираете каждый элемент массива. Вместо этого вы хотите перебирать только первые два измерения массива (ширина и высота).
Итерация по трехмерному массиву
В качестве быстрого примера, чтобы воспроизвести то, что вы видите без ArcMap:
import numpy as np
data = np.random.random((3, 10, 10))
for value in np.nditer(data):
print value
(Быстрое примечание: я использую arcpy форма съезда nbands x nrows x ncolumns Вот. Также очень часто можно увидеть nrows x ncolumns x nbands , В этом случае выражения индексации в последующих разделах будут другими)
Снова, nditer это не то, что вы хотите, поэтому, если вы действительно хотите сделать именно это (каждое значение в массиве вместо каждого r, g, b пикселя), это было бы гораздо более читабельным:
import numpy as np
data = np.random.random((3, 10, 10))
for value in data.flat:
print value
Эти два идентичны в этом случае.
Итерация по пикселям
Двигаясь дальше, вы хотите перебирать каждый пиксель. В этом случае вы бы сделали что-то вроде:
import numpy as np
data = np.random.random((3, 10, 10))
for pixel in data.reshape(3, -1).T:
r, g, b = pixel
print r, g, b
В этом случае мы временно рассматривали массив 10x10x3 как массив 100x3. Поскольку по умолчанию числовые массивы выполняют итерацию по первой оси, то итерация выполняется по каждому элементу r, g, b.
Если вы предпочитаете, вы также можете использовать индексирование напрямую, хотя это будет немного медленнее:
import numpy as np
data = np.random.random((3, 10, 10))
for i, j in np.ndindex(data.shape[:-2]):
r, g, b = data[:, i, j]
print r, g, b
Векторизация, не перебирать numpy массив
Тем не менее, в целом, повторение массива поэлементно, как это, не является эффективным способом использования numpy,
Вы упомянули, что пытаетесь определить, когда полосы были удалены и / или установлены на постоянное значение.
Вы можете иметь в виду три вещи: 1) есть только одна полоса, 2) данные в некоторых полосах были установлены на 0 (или другое значение), 3) изображение в оттенках серого, но хранится как RGB.
Вы можете проверить количество полос, посмотрев на массив numpy:
nbands = data.shape[0]
Или с помощью arcpy непосредственно:
nbands = raster.bandCount
Это касается первого случая, однако, похоже, что вы пытаетесь определить, когда у групп нет информации, в отличие от того, присутствуют они или нет.
Если вы всегда ожидаете иметь по крайней мере красный, зеленый и синий (иногда альфа, иногда нет), проще всего распаковать полосы, похожие на:
r, g, b = data[:3, :, :]
Таким образом, если есть альфа-группа, мы просто проигнорируем ее, а если ее нет, это не будет иметь значения. Опять же, это предполагает, что форма ваших данных - nbands x nrows x ncolumns (а не nrows x ncolumns x nbands).
Далее, если мы хотим проверить, все ли значения пикселей в полосе равны нулю, не перебираем. Вместо этого используйте numy логические сравнения. Они будут намного (>100x) быстрее:
r, g, b = data[:3, :, :]
print np.all(r == 0) # Are all red values zero?
Тем не менее, я думаю, что вы чаще всего хотите обнаружить это изображение в градациях серого, которое хранится как RGB. В этом случае значения красного, зеленого и синего каждого пикселя будут равны, но пиксели не будут одинаковыми. Вы можете проверить это, выполнив:
gray = (r == b) & (b == g)
print np.all(gray)
В общем, вы действительно не хотите перебирать каждый пиксель в массиве. Вместо этого используйте векторизованные выражения.
Предполагая, что вы уже знаете размер изображения (n x m), и ваш массив 1d numpy - A, это будет работать.
img2D = np.reshape(A, (m,n)).T
Пример: допустим, ваш массив изображений
img2D = array([[1, 2],
[3, 4],
[5, 6]])
Но вам дано A = array([1, 3, 5, 2, 4, 6]) Результат, который вы хотите получить:
img2D = np.reshape(A, (2, 3)).T