Fixing edge/pixel defects in image
What you are looking at for this problem is a combination of "Image Segmentation" and "Spatial Analysis".

The image segmentation part of this particular task is very simple because each pixel "class" is already labeled by a colour. You might need an additional labeling step because (for example) the blue class contains two disjoint areas, but this is easily dealt with via a simple connected component labeling. (We will come back to this in more detail).

The more challenging bit here is what corresponds to relationships expressed in natural language as "the blue areas inside the red areas", or "the black areas adjacent to red areas" and so on. These are relationships that correspond to mathematical "queries" such as "Point on line", "Point within area", "Area within area" and so on.

For example, to determine if a point is within a circle, all you have to test for is whether the distance from the centre of that circle to the point in question is less than the radius of the circle. To test a point within a rectangle you would test if its coordinates are within the boundaries of the rectangle and so on.

As you may have guessed, this does not sound easy to generalise (e.g. any kind of polygon, whether in 2D or 3D or higher dimensions etc). Fortuneatly, there are a lot of good libraries out there that implement geometrical primitives, operations and relationships for us to make such tasks easier.

In the case of Python there is Shapely. Shapely contains primitives such as Point, Line, Polygon and others as well as operators and relationships you can use to "query" their spatial relationships.

Here is a very simple example that answers part of your task: Is one area within another area?
Notice here, I simply wrote "area" as this is defined by an "unordered" set of pixels similar to what you would get out of simple queries such as "Give me all pixels in an image that are red".

    import shapely.geometry

    Q = shapely.geometry.MultiPoint([(u,v) for u in range(0,10) for v in range(0,10)])
    W = shapely.geometry.MultiPoint([(u,v) for u in range(5,8) for v in range(5,8)])


With these two lines we have created two MultiPoint (areas composed of many Points) "areas": One big one from (0,0) to (9,9) and one small one from (5,5) to (7,7). Furthermore, the small "patch" is contained by the bigger "patch". 

Let's probe that relationship

    Q.contains(W)


This should return True

   W.contains(Q)


This should return False.

You might be wondering what would happen if what we called the "small patch" contained two disjoint patches. In that case, Shapely would respond with False if only one of the patches was found within a testing area.

The black areas are a bit more challenging. You might think that their case is covered by the .touches() relationship but "touching" implies some sort of overlap and there are no criteria here by which to make those areas somehow to overlap. But what you can probe instead is distance.

In the case of MultiPoint objects, the distance between two objects is the shortest possible distance between their points. Again, here is a minimal example:

    import shapely.geometry

    Q = shapely.geometry.MultiPoint([(u,v) for u in range(0,10) for v in range(0,10)])
    W = shapely.geometry.MultiPoint([(u,v) for u in range(0,10) for v in range(12,22)])


Notice the differences in the way the "boxes" are setup here. To get the distance between them:

    Q.distance(W)


This should return 3.

You can see here that the black areas distance to blue, red and green is 1 (they are adjacent). Therefore, to locate which parts is a given black part adjacent to, you just assess its distance from the parts of interest.

So, with this bit out of the way, let's get back to how you would segment/label your pixel masks.

If your regions where not overlapping it would be extremely simple to extract the pixel coordinates that correspond to all pixels that have a specific color and create MultiPoint instances from them. Here is an example using Python's numpy and Shapely.

    import skimage.io
    import numpy
    import shapely.geometry

    # Read the image
    I = skimage.io.imread("someRRGBimage.png")
    # Exrtact the pixels that are equal to a specific value
    pixel_bunch = numpy.where((I[:,:,0]==the_red) &amp; (I[:,:,1]==the_green) &amp; (I[:,:,2]==the_blue)) # Notice here the_red, the_green, the_blue are the components of the area's colour.
    # Create a MultiPoint patch
    point_patch = shapely.geometry.MultiPoint(list(zip(pixel_bunch[0], pixel_bunch[1])));


But if you try to select the blues or the reds out of the image you are providing, you will get one large area, plus disjoint little ones and if you try to pass those through the spatial relationship analysis we are talking about above, they would fail.

The simplest way to split them apart is to use connectivity labeling. This works by scanning a binary image and assigning a different number to regions that are not connected. This basically inserts an extra step to the extraction of the pixel values above, as you would first have to create the mask image for a specific colour, then label that mask and finally create the regions based on the identified (and by now isolated) regions.

It's far easier to demonstrate this over grayscale images:

    import skimage.io
    import skimage.measure
    import numpy
    import Shapely.geometry

    # Read the image
    I = skimage.io.imread("someRRGBimage.png")
    # Work on the Red channel
    mask = (I[:,:,0]==255)
    labeled_mask = skimage.measure.label(mask)
    pixel_bunch = numpy.where(labeled_mask==some_patch_id) # Notice here some_patch_id is the zero based index of some isolated patch that has been labeled in the image.
    # Create a MultiPoint patch
    ...
    ...



It is not a difficult task but it sure is messy.

Finally, two notes:


You will probably have to re-colour those areas between the red and the blue that are currently black and leave black for the "background".
The complexity of deriving those relationships scales with the size of the regions, sometimes taking tens of seconds (or more) for very large MultiPoint sets. In that case, it might be more practical to work with polygons. This inserts yet another step just before deriving point_patch where you first derive the boundary of the identified pixel area and then instantiate a Polygon from that boundary. Operations between polygons will be much quicker.


Hope this helps.