Finding Raster Pixel Neighbors using GDAL?

I'm new to GDAL Python, and I'm working on a little recursive algorithm which requires finding the neighbors of a a pixel on a DEM raster. I've been able to do this by creating a little PixelLite class which is essentially just a tuple representing the raster coordinates of a pixel. Here's a bit of its implementation (feel free to skip):

``````    # a lighweight version of a pixel for quick computation
class PixelLite:

def __init__(self, x, y, raster_meta):
# set x and y
self.x = x
self.y = y

# set max_x and max_y
self.raster_meta = raster_meta
self.max_x = raster_meta.x_max
self.max_y = raster_meta.y_max

# set error flag
self.flag = False

# the set of all neighbors
self.neighbors = set()

if self.x > self.max_x | self.x < 0 | self.y > self.max_y | self.y < 0: # something has gone terribly wrong
self.flag = True
return

# !!! Below, I'm going to set out of bounds values to the value of the pixel itself, which will later be removed from the set !!!

# neighbor columns
if self.x > 0: # if left value will not be negative
left = self.x - 1 # add legit left
else:
left = self.x # add left value that will be removed

center = self.x

if self.x < self.max_x: # if right value will not be larger than raster dimensions
right = self.x + 1 # add legit right value
else:
right = self.x # add right value that will be removed

# Note: bottom and top are arbitrary names, what matters is that values do no exceed bounds; it may be the case that, because of raster indexing 'bottom' actually represents the pixel above the current pixel
# neighbor columns
if self.y > 0: # if bottom value will not be negative
bottom = self.y - 1 # add legit bottom value
else:
bottom = self.y # add bottom value that will be removed

middle = self.y

if self.y < self.max_y: # if top value will not exceed raster extent
top = self.y + 1 # add legit top value
else:
top = self.y # add top value that will be removed

# top row of neighbors

# middle row of neighbors

# bottom row of neighbors

# the method I'm using will sometimes add the center, middle element as I've used that as a default to set values to if there is an out of bounds error. So let's just pop that off...; a pixel cannot be its own neighbor
Is there a faster `getPixelNeighbors` function provided by GDAL or another raster library i.e. one which would allow me to do the same but which might have some neat little optimization tricks up its sleeve?
• There may be optimised neighbourhood code buried in the GDAL C++ lbrary for use by other function, but I'm not aware of it. I would use the gdal `Dataset` or `RasterBand` `ReadAsArray` to read as a numpy array and use the power of numpy indexing. Even easier, use rasterio, your open dataset object has a `read` method for reading as a numpy array and an `index` method to transform map coordinates to pixel indexes which you can use with your numpy array. – user2856 Apr 20 at 3:05