Merge pull request #712 from wiredfool/morphology

Morphology

PIL/ImageMorph.py

@@ -0,0 +1,236 @@
+# A binary morphology add-on for the Python Imaging Library
+#
+# History:
+#   2014-06-04 Initial version.
+#
+# Copyright (c) 2014 Dov Grobgeld <dov.grobgeld@gmail.com>
+
+from PIL import Image
+from PIL import _imagingmorph
+import re
+
+LUT_SIZE = 1<<9
+class LutBuilder:
+    """A class for building MorphLut's from a descriptive language
+
+      The input patterns is a list of a strings sequences like these:
+
+      4:(...
+         .1.
+         111)->1
+
+      (whitespaces including linebreaks are ignored). The option 4
+      descibes a series of symmetry operations (in this case a
+      4-rotation), the pattern is decribed by:
+
+         . or X - Ignore
+         1 - Pixel is on
+         0 - Pixel is off
+
+      The result of the operation is described after "->" string.
+
+      The default is to return the current pixel value, which is
+      returned if no other match is found.
+
+      Operations:
+         4 - 4 way rotation
+         N - Negate
+         1 - Dummy op for no other operation (an op must always be given)
+         M - Mirroring
+
+      Example:
+
+      lb = LutBuilder(patterns = ["4:(... .1. 111)->1"])
+      lut = lb.build_lut()
+         
+    """
+    def __init__(self,patterns = None,op_name=None):
+        if patterns is not None:
+            self.patterns = patterns
+        else:
+            self.patterns = []
+        self.lut = None
+        if op_name is not None:
+            known_patterns = {
+                'corner' : ['1:(... ... ...)->0',
+                            '4:(00. 01. ...)->1'],
+                'dilation4' : ['4:(... .0. .1.)->1'],
+                'dilation8' : ['4:(... .0. .1.)->1',
+                               '4:(... .0. ..1)->1'],
+                'erosion4' : ['4:(... .1. .0.)->0'],
+                'erosion8' : ['4:(... .1. .0.)->0',
+                              '4:(... .1. ..0)->0'],
+                'edge' : ['1:(... ... ...)->0',
+                          '4:(.0. .1. ...)->1',
+                          '4:(01. .1. ...)->1']
+            }
+            if not op_name in known_patterns:
+                raise Exception('Unknown pattern '+op_name+'!')
+
+            self.patterns = known_patterns[op_name]
+
+    def add_patterns(self, patterns):
+        self.patterns += patterns
+
+    def build_default_lut(self):
+        symbols = [0, 1]
+        m = 1 << 4  # pos of current pixel
+        self.lut = bytearray([symbols[(i & m)>0] for i in range(LUT_SIZE)])
+
+    def get_lut(self):
+        return self.lut
+
+    def _string_permute(self, pattern, permutation):
+        """string_permute takes a pattern and a permutation and returns the
+        string permuted accordinging to the permutation list.
+        """
+        assert(len(permutation)==9)
+        return ''.join([pattern[p] for p in permutation])
+
+    def _pattern_permute(self, basic_pattern, options, basic_result):
+        """pattern_permute takes a basic pattern and its result and clones
+        the mattern according to the modifications described in the $options
+        parameter. It returns a list of all cloned patterns."""
+        patterns = [(basic_pattern, basic_result)]
+
+        # rotations
+        if '4' in options:
+            res = patterns[-1][1]
+            for i in range(4):
+                patterns.append(
+                    (self._string_permute(patterns[-1][0],
+                                    [6,3,0,
+                                     7,4,1,
+                                     8,5,2]), res))
+        # mirror
+        if 'M' in options:
+            n = len(patterns)
+            for pattern,res in patterns[0:n]:
+                patterns.append(
+                    (self._string_permute(pattern, [2,1,0,
+					      5,4,3,
+					      8,7,6]), res))
+
+        # negate
+        if 'N' in options:
+            n = len(patterns)
+            for pattern,res in patterns[0:n]:
+                # Swap 0 and 1
+                pattern = (pattern
+                           .replace('0','Z')
+                           .replace('1','0')
+                           .replace('Z','1'))
+                res = '%d'%(1-int(res))
+                patterns.append((pattern, res))
+
+        return patterns
+
+    def build_lut(self):
+        """Compile all patterns into a morphology lut.
+
+        TBD :Build based on (file) morphlut:modify_lut
+        """
+        self.build_default_lut()
+        patterns = []
+
+        # Parse and create symmetries of the patterns strings
+        for p in self.patterns:
+            m = re.search(r'(\w*):?\s*\((.+?)\)\s*->\s*(\d)', p.replace('\n',''))
+            if not m:
+                raise Exception('Syntax error in pattern "'+p+'"')
+            options = m.group(1)
+            pattern = m.group(2)
+            result = int(m.group(3))
+
+            # Get rid of spaces
+            pattern= pattern.replace(' ','').replace('\n','')
+
+            patterns += self._pattern_permute(pattern, options, result)
+
+#        # Debugging
+#        for p,r in patterns:
+#            print p,r
+#        print '--'
+
+        # compile the patterns into regular expressions for speed
+        for i in range(len(patterns)):
+            p = patterns[i][0].replace('.','X').replace('X','[01]')
+            p = re.compile(p)
+            patterns[i] = (p, patterns[i][1])
+
+        # Step through table and find patterns that match.
+        # Note that all the patterns are searched. The last one
+        # caught overrides
+        for i in range(LUT_SIZE):
+            # Build the bit pattern
+            bitpattern = bin(i)[2:]
+            bitpattern = ('0'*(9-len(bitpattern)) + bitpattern)[::-1]
+
+            for p,r in patterns:
+                if p.match(bitpattern):
+                    self.lut[i] = [0, 1][r]
+
+        return self.lut
+        
+class MorphOp:
+    """A class for binary morphological operators"""
+
+    def __init__(self,
+                 lut=None,
+                 op_name = None,
+                 patterns = None):
+        """Create a binary morphological operator"""
+        self.lut = lut
+        if op_name is not None:
+            self.lut = LutBuilder(op_name = op_name).build_lut()
+        elif patterns is not None:
+            self.lut = LutBuilder(patterns = patterns).build_lut()
+
+    def apply(self, image):
+        """Run a single morphological operation on an image
+
+        Returns a tuple of the number of changed pixels and the
+        morphed image"""
+        if self.lut is None:
+            raise Exception('No operator loaded')
+        
+        outimage = Image.new(image.mode, image.size, None)
+        count = _imagingmorph.apply(bytes(self.lut), image.im.id, outimage.im.id)
+        return count, outimage
+    
+    def match(self, image):
+        """Get a list of coordinates matching the morphological operation on an image
+
+        Returns a list of tuples of (x,y) coordinates of all matching pixels."""
+        if self.lut is None:
+            raise Exception('No operator loaded')
+        
+        return  _imagingmorph.match(bytes(self.lut), image.im.id)
+
+    def get_on_pixels(self, image):
+        """Get a list of all turned on pixels in a binary image
+        Returns a list of tuples of (x,y) coordinates of all matching pixels."""
+        
+        return  _imagingmorph.get_on_pixels(image.im.id)
+
+    def load_lut(self, filename):
+        """Load an operator from an mrl file"""
+        with open(filename,'rb') as f:
+            self.lut = bytearray(f.read())
+            
+        if len(self.lut)!= 8192:
+            self.lut = None
+            raise Exception('Wrong size operator file!')
+
+    def save_lut(self, filename):
+        """Load an operator save mrl file"""
+        if self.lut is None:
+            raise Exception('No operator loaded')
+        with open(filename,'wb') as f:
+            f.write(self.lut)
+
+    def set_lut(self, lut):
+        """Set the lut from an external source"""
+        self.lut = lut
+
+    

Tests/test_imagemorph.py

@@ -0,0 +1,169 @@
+# Test the ImageMorphology functionality
+from helper import *
+
+from PIL import Image
+from PIL import ImageMorph
+
+
+class MorphTests(PillowTestCase):
+
+    def setUp(self):
+        self.A = self.string_to_img(
+            """    
+            .......
+            .......
+            ..111..
+            ..111..
+            ..111..
+            .......
+            .......
+            """
+            )    
+    
+
+    def img_to_string(self, im):
+        """Turn a (small) binary image into a string representation"""
+        chars = '.1'
+        width, height = im.size
+        return '\n'.join(
+            [''.join([chars[im.getpixel((c,r))>0] for c in range(width)])
+             for r in range(height)])
+
+    def string_to_img(self, image_string):
+        """Turn a string image representation into a binary image"""
+        rows = [s for s in image_string.replace(' ','').split('\n')
+                if len(s)]
+        height = len(rows)
+        width = len(rows[0])
+        im = Image.new('L',(width,height))
+        for i in range(width):
+            for j in range(height):
+                c = rows[j][i]
+                v = c in 'X1'
+                im.putpixel((i,j),v)
+
+        return im
+
+    def img_string_normalize(self, im):
+        return self.img_to_string(self.string_to_img(im))
+
+    def assert_img_equal(self, A, B):
+        self.assertEqual(self.img_to_string(A), self.img_to_string(B))
+
+    def assert_img_equal_img_string(self, A, Bstring):
+        self.assertEqual(self.img_to_string(A), self.img_string_normalize(Bstring))
+
+    def test_str_to_img(self):
+        im = Image.open('Tests/images/morph_a.png')
+        self.assert_image_equal(self.A, im)
+
+    def create_lut(self):
+        for op in ('corner', 'dilation4', 'dilation8', 'erosion4', 'erosion8', 'edge'):
+            lb = ImageMorph.LutBuilder(op_name=op)
+            lut = lb.build_lut(self)
+            with open('Tests/images/%s.lut' % op, 'wb') as f:
+                f.write(lut)
+
+    #create_lut()
+    def test_lut(self):
+        for op in ('corner', 'dilation4', 'dilation8', 'erosion4', 'erosion8', 'edge'):
+            lb = ImageMorph.LutBuilder(op_name=op)
+            lut = lb.build_lut()
+            with open('Tests/images/%s.lut' % op , 'rb') as f:
+                self.assertEqual(lut, bytearray(f.read()))
+
+
+    # Test the named patterns
+    def test_erosion8(self):
+        # erosion8
+        mop = ImageMorph.MorphOp(op_name='erosion8')
+        count,Aout = mop.apply(self.A)
+        self.assertEqual(count,8)
+        self.assert_img_equal_img_string(Aout,
+                                         """
+                                         .......
+                                         .......
+                                         .......
+                                         ...1...
+                                         .......
+                                         .......
+                                         .......
+                                         """)
+
+    def test_dialation8(self):
+        # dialation8
+        mop = ImageMorph.MorphOp(op_name='dilation8')
+        count,Aout = mop.apply(self.A)
+        self.assertEqual(count,16)
+        self.assert_img_equal_img_string(Aout,
+                                         """
+                                         .......
+                                         .11111.
+                                         .11111.
+                                         .11111.
+                                         .11111.
+                                         .11111.
+                                         .......
+                                         """)
+
+    def test_erosion4(self):
+        # erosion4
+        mop = ImageMorph.MorphOp(op_name='dilation4')
+        count,Aout = mop.apply(self.A)
+        self.assertEqual(count,12)
+        self.assert_img_equal_img_string(Aout,
+                                         """
+                                         .......
+                                         ..111..
+                                         .11111.
+                                         .11111.
+                                         .11111.
+                                         ..111..
+                                         .......
+                                         """)
+
+    def test_edge(self):
+        # edge 
+        mop = ImageMorph.MorphOp(op_name='edge')
+        count,Aout = mop.apply(self.A)
+        self.assertEqual(count,1)
+        self.assert_img_equal_img_string(Aout,
+                                         """
+                                         .......
+                                         .......
+                                         ..111..
+                                         ..1.1..
+                                         ..111..
+                                         .......
+                                         .......
+                                         """)
+
+
+    def test_corner(self):
+        # Create a corner detector pattern
+        mop = ImageMorph.MorphOp(patterns = ['1:(... ... ...)->0',
+                                             '4:(00. 01. ...)->1'])
+        count,Aout = mop.apply(self.A)
+        self.assertEqual(count,5)
+        self.assert_img_equal_img_string(Aout,
+                                         """
+                                         .......
+                                         .......
+                                         ..1.1..
+                                         .......
+                                         ..1.1..
+                                         .......
+                                         .......
+                                         """)
+        
+
+        # Test the coordinate counting with the same operator
+        coords = mop.match(self.A)
+        self.assertEqual(len(coords), 4)
+        self.assertEqual(tuple(coords),
+                     ((2,2),(4,2),(2,4),(4,4)))
+
+        coords = mop.get_on_pixels(Aout)
+        self.assertEqual(len(coords), 4)
+        self.assertEqual(tuple(coords),
+                     ((2,2),(4,2),(2,4),(4,4)))

_imagingmorph.c

@@ -0,0 +1,303 @@
+/*
+ * The Python Imaging Library
+ *
+ * A binary morphology add-on for the Python Imaging Library
+ *
+ * History:
+ *   2014-06-04 Initial version.
+ *
+ * Copyright (c) 2014 Dov Grobgeld <dov.grobgeld@gmail.com>
+ *
+ * See the README file for information on usage and redistribution.
+ */
+
+#include "Python.h"
+#include "Imaging.h"
+#include "py3.h"
+
+#define LUT_SIZE (1<<9)
+
+/* Apply a morphologic LUT to a binary image. Outputs a
+   a new binary image.
+
+   Expected parameters:
+
+      1. a LUT - a 512 byte size lookup table.
+      2. an input Imaging image id.
+      3. an output Imaging image id
+
+   Returns number of changed pixels.
+*/
+static PyObject*
+apply(PyObject *self, PyObject* args)
+{
+    const char *lut;
+    PyObject *py_lut;
+    Py_ssize_t lut_len, i0, i1;
+    Imaging imgin, imgout;
+    int width, height;
+    int row_idx, col_idx;
+    UINT8 **inrows, **outrows;
+    int num_changed_pixels = 0;
+
+    if (!PyArg_ParseTuple(args, "Onn", &py_lut, &i0, &i1)) {
+        PyErr_SetString(PyExc_RuntimeError, "Argument parsing problem");
+        return NULL;
+    }
+
+    if (!PyBytes_Check(py_lut)) {
+        PyErr_SetString(PyExc_RuntimeError, "The morphology LUT is not a bytes object");
+        return NULL;
+    }
+
+    lut_len = PyBytes_Size(py_lut);
+
+    if (lut_len < LUT_SIZE) {
+        PyErr_SetString(PyExc_RuntimeError, "The morphology LUT has the wrong size");
+        return NULL;
+    }
+
+    lut = PyBytes_AsString(py_lut);
+
+    imgin = (Imaging) i0;
+    imgout = (Imaging) i1;
+    width = imgin->xsize;
+    height = imgin->ysize;
+
+    if (imgin->type != IMAGING_TYPE_UINT8 &&
+        imgin->bands != 1) {
+        PyErr_SetString(PyExc_RuntimeError, "Unsupported image type");
+        return NULL;
+    }
+    if (imgout->type != IMAGING_TYPE_UINT8 &&
+        imgout->bands != 1) {
+        PyErr_SetString(PyExc_RuntimeError, "Unsupported image type");
+        return NULL;
+    }
+
+    inrows = imgin->image8;
+    outrows = imgout->image8;
+
+    for (row_idx=0; row_idx < height; row_idx++) {
+        UINT8 *outrow = outrows[row_idx];
+        UINT8 *inrow = inrows[row_idx];
+        UINT8 *prow, *nrow; /* Previous and next row */
+
+        /* zero boundary conditions. TBD support other modes */
+        outrow[0] = outrow[width-1] = 0;
+        if (row_idx==0 || row_idx == height-1) {
+            for(col_idx=0; col_idx<width; col_idx++)
+                outrow[col_idx] = 0;
+            continue;
+        }
+
+        prow = inrows[row_idx-1];
+        nrow = inrows[row_idx+1];
+
+        for (col_idx=1; col_idx<width-1; col_idx++) {
+            int cim = col_idx-1;
+            int cip = col_idx+1;
+            unsigned char b0 = prow[cim] &1;
+            unsigned char b1 = prow[col_idx]&1;
+            unsigned char b2 = prow[cip]&1;
+
+            unsigned char b3 = inrow[cim]&1;
+            unsigned char b4 = inrow[col_idx]&1;
+            unsigned char b5 = inrow[cip]&1;
+
+            unsigned char b6 = nrow[cim]&1;
+            unsigned char b7 = nrow[col_idx]&1;
+            unsigned char b8 = nrow[cip]&1;
+
+            int lut_idx = (b0 
+                           |(b1 << 1)
+                           |(b2 << 2)
+                           |(b3 << 3)
+                           |(b4 << 4)
+                           |(b5 << 5)
+                           |(b6 << 6)
+                           |(b7 << 7)
+                           |(b8 << 8));
+            outrow[col_idx] = 255*(lut[lut_idx]&1);
+            num_changed_pixels += ((b4&1)!=(outrow[col_idx]&1));
+        }
+    }
+    return Py_BuildValue("i",num_changed_pixels);
+}
+
+/* Match a morphologic LUT to a binary image and return a list
+   of the coordinates of all all matching pixels.
+
+   Expected parameters:
+
+      1. a LUT - a 512 byte size lookup table.
+      2. an input Imaging image id.
+
+   Returns list of matching pixels.
+*/
+static PyObject*
+match(PyObject *self, PyObject* args)
+{
+    const char *lut;
+    PyObject *py_lut;
+    Py_ssize_t lut_len, i0;
+    Imaging imgin;
+    int width, height;
+    int row_idx, col_idx;
+    UINT8 **inrows;
+    PyObject *ret = PyList_New(0);
+
+    if (!PyArg_ParseTuple(args, "On", &py_lut, &i0)) {
+        PyErr_SetString(PyExc_RuntimeError, "Argument parsing problem");
+        return NULL;
+    }
+
+    if (!PyBytes_Check(py_lut)) {
+        PyErr_SetString(PyExc_RuntimeError, "The morphology LUT is not a bytes object");
+        return NULL;
+    }
+
+    lut_len = PyBytes_Size(py_lut);
+
+    if (lut_len < LUT_SIZE) {
+        PyErr_SetString(PyExc_RuntimeError, "The morphology LUT has the wrong size");
+        return NULL;
+    }
+
+    lut = PyBytes_AsString(py_lut);
+    imgin = (Imaging) i0;
+
+    if (imgin->type != IMAGING_TYPE_UINT8 &&
+        imgin->bands != 1) {
+        PyErr_SetString(PyExc_RuntimeError, "Unsupported image type");
+        return NULL;
+    }
+
+    inrows = imgin->image8;
+    width = imgin->xsize;
+    height = imgin->ysize;
+
+    for (row_idx=1; row_idx < height-1; row_idx++) {
+        UINT8 *inrow = inrows[row_idx];
+        UINT8 *prow, *nrow;
+
+        prow = inrows[row_idx-1];
+        nrow = inrows[row_idx+1];
+
+        for (col_idx=1; col_idx<width-1; col_idx++) {
+            int cim = col_idx-1;
+            int cip = col_idx+1;
+            unsigned char b0 = prow[cim] &1;
+            unsigned char b1 = prow[col_idx]&1;
+            unsigned char b2 = prow[cip]&1;
+
+            unsigned char b3 = inrow[cim]&1;
+            unsigned char b4 = inrow[col_idx]&1;
+            unsigned char b5 = inrow[cip]&1;
+
+            unsigned char b6 = nrow[cim]&1;
+            unsigned char b7 = nrow[col_idx]&1;
+            unsigned char b8 = nrow[cip]&1;
+
+            int lut_idx = (b0 
+                           |(b1 << 1)
+                           |(b2 << 2)
+                           |(b3 << 3)
+                           |(b4 << 4)
+                           |(b5 << 5)
+                           |(b6 << 6)
+                           |(b7 << 7)
+                           |(b8 << 8));
+            if (lut[lut_idx]) {
+                PyObject *coordObj = Py_BuildValue("(nn)",col_idx,row_idx);
+                PyList_Append(ret, coordObj);
+            }
+        }
+    }
+
+    return ret;
+}
+
+/* Return a list of the coordinates of all turned on pixels in an image.
+   May be used to extract features after a sequence of MorphOp's were applied.
+   This is faster than match as only 1x1 lookup is made.
+*/
+static PyObject*
+get_on_pixels(PyObject *self, PyObject* args)
+{
+    Py_ssize_t i0;
+    Imaging img;
+    UINT8 **rows;
+    int row_idx, col_idx;
+    int width, height;
+    PyObject *ret = PyList_New(0);
+  
+    if (!PyArg_ParseTuple(args, "n", &i0)) {
+        PyErr_SetString(PyExc_RuntimeError, "Argument parsing problem");
+
+        return NULL;
+    }
+    img = (Imaging) i0;
+    rows = img->image8;
+    width = img->xsize;
+    height = img->ysize;
+  
+    for (row_idx=0; row_idx < height; row_idx++) {
+        UINT8 *row = rows[row_idx];
+        for (col_idx=0; col_idx<width; col_idx++) {
+            if (row[col_idx]) {
+                PyObject *coordObj = Py_BuildValue("(nn)",col_idx,row_idx);
+                PyList_Append(ret, coordObj);
+            }
+        }
+    }
+    return ret;
+}
+
+
+static int
+setup_module(PyObject* m)
+{
+    PyObject* d = PyModule_GetDict(m);
+
+    PyDict_SetItemString(d, "__version", PyUnicode_FromString("0.1"));
+
+    return 0;
+}
+
+static PyMethodDef functions[] = {
+    /* Functions */
+    {"apply", (PyCFunction)apply, 1},
+    {"get_on_pixels", (PyCFunction)get_on_pixels, 1},
+    {"match", (PyCFunction)match, 1},
+};
+
+#if PY_VERSION_HEX >= 0x03000000
+PyMODINIT_FUNC
+PyInit__imagingmorph(void) {
+    PyObject* m;
+
+    static PyModuleDef module_def = {
+        PyModuleDef_HEAD_INIT,
+        "_imagingmorph",         /* m_name */
+        "A module for doing image morphology",               /* m_doc */
+        -1,                 /* m_size */
+        functions,          /* m_methods */
+    };
+
+    m = PyModule_Create(&module_def);
+
+    if (setup_module(m) < 0)
+        return NULL;
+
+    return m;
+}
+#else
+PyMODINIT_FUNC
+init_imagingmorph(void)
+{
+    PyObject* m = Py_InitModule("_imagingmorph", functions);
+    setup_module(m);
+}
+#endif
+

setup.py

@@ -571,6 +571,9 @@ class pil_build_ext(build_ext):
         if os.path.isfile("_imagingmath.c"):
             exts.append(Extension("PIL._imagingmath", ["_imagingmath.c"]))
 
+        if os.path.isfile("_imagingmorph.c"):
+            exts.append(Extension("PIL._imagingmorph", ["_imagingmorph.c"]))
+
         self.extensions[:] = exts
 
         build_ext.build_extensions(self)