Vjylku

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

            self.neighbors.add((left, top))

            self.neighbors.add((center, top))

            self.neighbors.add((right, top))



            # middle row of neighbors

            self.neighbors.add((left, middle))

            self.neighbors.add((right, middle))



            # bottom row of neighbors

            self.neighbors.add((left, bottom))

            self.neighbors.add((center, bottom))

            self.neighbors.add((right, bottom))



            # 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

            self.neighbors.discard((center, middle))



        ...

Basically, what this does is create the tuples representing each neighbor for a given pixel whenever the PixelLite class is instantiated. Although this is pretty lightweight, because it is implemented in very basic Python without any kind of binding, I would imagine it lacks the optimization of a more sophisticated solution such as I would expect to find in GDAL.

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?