FINDING OPTIMAL SEQUENTIAL DECOMPOSITIONS OF EROSIONS AND DILATIONS

A. JUNIOR BARRERA, B. CARLOS E. FERREIRA and

C. RONALDO F. HASHIMOTO ^{∗}

Departamento de Ciência da Computação do IME/USP Rua do Matão, 1010, CEP 05508-900, São Paulo, SP, Brasil. email: [jb, cef, ronaldo]@ime.usp.br

Abstract. A fundamental problem in the field of automatic design of morphological operators is the transformation of standard morphological representations, that are adequate for operator learning, into other equivalent representations, that permit efficient implementations of the operators learned. In this paper, we study this problem restricted to the families of set translation invariant erosions and dilations. We present a general method for the decomposition of erosions (respec., dilations) by arbitrarily shaped structuring elements into a minimum number of compositions of erosions (respec. dilations) by subsets of the elementary square (i.e.,the 3 × 3 square centered at the origin). The fundamental idea of the method proposed is the application of Combinatorial Optimization techniques, supported by generic algebraic and geometrical properties of erosions and dilations.

Key words: erosion, dilation, Minkowski addition, decomposition, structuring element.

1. Introduction

The problem of designing automatically morphological operators [1] may be approached by a process composed of two main steps: i - learning of the target operator by the estimation of the structuring elements (SE’s) that characterize it in a standard representation form; ii - transformation of the standard form morphological operator learned into other equivalent representations that permit more efficient implementations.

We find in the literature some works that study the problem of transformation of decomposition structures [3, 2]. However, the general problem of transforming the sup or inf-representation [3] of t.i. set operators into efficient sequential or hybrid (sequential-parallel) representations is extremely hard and practically not studied. In fact, some aspects of this problem have been studied for the families of erosions and dilations.

The speed up achieved by representing erosions and dilations by sequential decompositions, in conventional machines, was quantitatively studied by Maragos [5, p. 77], who showed examples where the complexity of the algorithms that implement erosions and dilations went from quadratic, in the usual representation (i.e., the sup or inf-representation), to linear, in the sequential decomposition.

∗ The authors are presented in alphabetic order.

In fact, the problem of transforming the usual representation of an erosion or a dilation into a sequential representation is equivalent to find a decomposition for the SE that characterizes it [9, p. 47].

Several researchers [8, 10, 6] have studied the problem of decomposing SE’s as sequences of Minkowski additions of smaller subsets and proposed different algorithms to generate them. Xu [10] developed an algorithm for the decomposition of convex SE’s in terms of subsets of the elementary square (i.e., the 3 × 3 square centered at the origin). Park and Chin [6] developed an extension of Xu’s algorithm for the decomposition of simply connected (i.e., an 8-connected set that contains no holes) SE’s. Both algorithms are optimal, in the sense that they produce decompositions that use the minimum number of subsets.

We should remark that not all SE’s have sequential decompositions by Minkowski additions [10]. Furthermore, it is not known an efficient algorithm for determining the existence of such decompositions for an arbitrary SE.

Here, we present a method for the generation of decompositions of arbitrarily shaped SE’s as sequences of Minkowski additions of subsets of the elementary square. If there exist such decompositions, the method gives the optimal one, otherwise, it proves that the SE is not decomposable. The formulation adopted is a branch and bound search in a tree that represents the space of all possible sequences of subsets.

Following this introduction, Section 2 presents definitions and properties of Minkowski additions. Section 3 recalls the formulation of Combinatorial Optimization problems by the branch and bound technique. Section 4 joins the information of Sections 2 and 3 for presenting the method proposed. Section 5 presents some experimental results. Finally, Section 6 discusses the main contributions of this paper and some possible next steps of this research.

2. Definitions and Properties

For any X ⊂ ZZ^{2} and y ∈ ZZ^{2}, X_{y} denotes the translation of X by y, that is, X_{y} = {x ∈ ZZ^{2} : x − y ∈ X}.

Let X, Y ⊂ ZZ^{2}. The Minkowski addition and subtraction of X and Y are the subsets given, respectively, by X⊕Y = ∪{X_{y} : y ∈ Y } and X⊖Y = ∩{X_{−y} : y ∈ Y }.

Two properties of Minkowski addition are commutativity (i.e., X ⊕ Y = Y ⊕ X) and translation by unitary sets (i.e., X ⊕ {y} = X_{y}).

In this paper, we will say that a SE A is decomposable by Minkowski additions if there exists a sequence of subsets B_{1}, …, B_{n} of the elementary square such that A = B_{1} ⊕ B_{2} ⊕ ⋯ ⊕ B_{n}.

The sequence of subsets that represent a SE can be decomposed in two subsequences: shape and translation. The shape sequence represents the shape of the SE, while the translation sequence defines its position in the integer plane. The shape sequence B_{1}, …, B_{k} and the number k are called, respectively, a decomposition of the SE A and the length of this decomposition.

Let A, B ⊂ ZZ^{2}. We say that B is an invariant of A if and only if (iff) A =

(A⊖B)⊕B. For example, the subset [ 0 1 1 ] is an invariant of the subset [ \begin{matrix}1 & 1 & 0 \\ 0 & 1 & 1 \end{matrix} ],

where the 1’s represent the points that belong to the subset, the 0’s represent points that do not belong to the subset and the bold face character represents the origin (i.e., the point (0, 0)).

Proposition 1 Let A, X, Y ⊂ ZZ^{2}. If A = X ⊕Y, then X and Y are both invariants of A.

In the following, we give two immediate corollaries of this proposition.

Corollary 2 If A = B_{1} ⊕ B_{2} ⊕ ⋯ ⊕ B_{k}, then each B_{i} is an invariant of A.

This corollary implies that every subset in a decomposition of a SE is an invariant of it. So, to find a decomposition of a SE, we can restrict our attention only to its invariants.

Corollary 3 Let X, A ⊂ ZZ^{2}. If X is not an invariant of A, then, for any non empty subset Y ⊂ ZZ^{2}, Y ⊕ X is not an invariant of A.

Let B ⊂ ZZ^{2} and n be a positive integer. The succession of n − 1 self Minkowski additions of B is the subset denoted nB, that is, nB = ((B ⊕ B) ⊕ ⋯ ⊕ B).

Let A, B ⊂ ZZ^{2} such that B is an invariant of A. The multiplicity of B is the greatest positive integer n such that nB is an invariant of A. For example, the

multiplicity of the invariant B = [ 0 1 1 ] of A = [ \begin{matrix}1 & 1 & 1 \\ 1 & 1 & 1 \end{matrix} ] is 2, since 2B = [ 1 1 1 ]

is an invariant of A and 3B = [ 1 1 1 1 ] is not. Note that unitary sets have infinity multiplicity.

Proposition 4 Let A, B ⊂ ZZ^{2} such that |B| ≥ 2 and B is an invariant of A. If n is the multiplicity of B, then any decomposition of A contains at most n elements equal to B.

Let A ⊂ ZZ^{2} and B(A) = {B ⊆ {−1, 0, 1}^{2} : B is an invariant of A and |B| ≥ 2}. Let B_{1}, B_{2}, ⋯, B_{k} be all elements of B(A) in a fixed order and n_{i} be the multiplicity of the invariant B_{i}, i = 1, …, k. The invariant sequence of A is the sequence

\underbrace{B_{1}, ⋯, B_{1}}_{n_{1} times}, ⋯, \underbrace{B_{k}, ⋯, B_{k}}_{n_{k} times}.

Now, we will study the length of decompositions of an arbitrarily shaped SE from the analysis of measures taken over it.

Let u⃗_{1}, u⃗_{2}, u⃗_{3} and u⃗_{4} be the Cartesian axis that intersect the origin and have slopes, respectively, 0, 45, 90 and 135 degrees. For a given point x ∈ ZZ^{2}, let l_{1}(x), l_{2}(x), l_{3}(x), l_{4}(x) be the orthogonal projections of x at the Cartesian axis u⃗_{1}, u⃗_{2}, u⃗_{3} and u⃗_{4}, respectively. Observe that these projections are integer numbers at the directions 0 and 90 degrees and real numbers (proportional to \frac{\sqrt{2}}{2}) at the directions 45 and 135 degrees. Let x_{1}, x_{2}, x_{3} and x_{4} denote the normalized orthogonal projections of the point x at the Cartesian axis u⃗_{1}, u⃗_{2}, u⃗_{3} and u⃗_{4}, respectively, that is, x_{1} = l_{1}(x), x_{2} = l_{2}(x) · \sqrt{2}, x_{3} = l_{3}(x) and x_{4} = l_{4}(x) · \sqrt{2}. For example, the normalized projections of the point x = (−5, 2) ∈ ZZ^{2} are: x_{1} = −5, x_{2} = −3, x_{3} = 2 and x_{4} = 7.

Let A ⊂ ZZ^{2}. For i = 1, 2, 3, 4, let x_{i}^{max}(A) and x_{i}^{min}(A) be, respectively, the maximum and the minimum normalized orthogonal projection at the Cartesian axis u⃗_{i} of the points in A, that is, x_{i}^{max}(A) = max{x_{i} : x ∈ A} and x_{i}^{min}(A) = min{x_{i} : x ∈ A}. The orthogonal projection of the set A is the vector ρ(A) ∈ ZZ^{4} such that, for i = 1, 2, 3, 4, ρ_{i}(A) = x_{i}^{max}(A) − x_{i}^{min}(A). For example, the orthogonal projection

of the subset A = [ \begin{matrix}1 & 1 & 0 \\ 0 & 1 & 1 \end{matrix} ] is ρ(A) = (2, 1, 1, 3).

Proposition 5 If A ⊂ ZZ^{2}, then max{ρ_{2}(A)/2, ρ_{4}(A)/2} ≤ max{ρ_{1}(A), ρ_{3}(A)}.

Proof. ρ_{1}(A) and ρ_{3}(A) are the length of the sides of the smallest rectangle that contains A. It is easy to see that ρ_{2}(A)\frac{\sqrt{2}}{2} is less or equal to the length of a diagonal of this rectangle. Thus, ρ_{2}(A)^{2}/2 ≤ ρ_{1}(A)^{2} + ρ_{3}(A)^{2} ≤ 2(max{ρ_{1}(A), ρ_{3}(A)})^{2}. Hence, ρ_{2}(A)/2 ≤ max{ρ_{1}(A), ρ_{3}(A)}. Similarly, ρ_{4}(A)/2 ≤ max{ρ_{1}(A), ρ_{3}(A)}. Therefore, we have that max{ρ_{2}(A)/2, ρ_{4}(A)/2} ≤ max{ρ_{1}(A), ρ_{3}(A)}.

The following proposition gives a lower bound for the length of decompositions of a SE by Minkowski additions.

Proposition 6 If A ⊂ ZZ^{2} is decomposable, then a decomposition of A contains at least ⌈max{ρ_{1}(A), ρ_{3}(A)}/2⌉ subsets of the elementary square.

Proof. Let m be the number of subsets of the elementary square of a decomposition of the SE A. Clearly, A ⊆ mS, where S is the elementary square, and the orthogonal projection of mS is ρ(mS) = (2m, 4m, 2m, 4m). Hence, since A ⊆ mS, we have that ρ_{1}(A) ≤ 2m, ρ_{2}(A) ≤ 4m, ρ_{3}(A) ≤ 2m and ρ_{4}(A) ≤ 4m. Thus, max{ρ_{1}(A), ρ_{2}(A)/2, ρ_{3}(A), ρ_{4}(A)/2} ≤ 2m and, by Proposition 5, we have that max{ρ_{1}(A), ρ_{3}(A)} ≤ 2m. Therefore, m ≥ ⌈max{ρ_{1}(A), ρ_{3}(A)}/2⌉.

A sequence of subsets of the elementary square B_{1}, B_{2}, …, B_{k} is said an unfeasible sequence for the decomposition of a SE A iff the subset B_{1} ⊕ B_{2} ⊕ ⋯ ⊕ B_{k} is not an invariant of A. The next corollary is an immediate consequence of Proposition 1.

Corollary 7 If B_{1}, B_{2}, …, B_{k} is an unfeasible sequence of a given SE A, then B_{1}, B_{2}, …, B_{k} is not a subsequence of a decomposition of A.

In the following, we state some conditions that characterize unfeasible sequences for the decomposition of SE’s.

The convex hull C(A) of a set A ⊂ ZZ^{2} is the intersection of all the half planes that contains A. We will suppose that the subsets are represented in the square grid and consider just the half planes with slopes 0, 45, 90 and 135 degrees to build the convex hull. A set X ⊂ ZZ^{2} is said convex iff X = C(X).

Now, we define the projection vector of a SE A as a point in ZZ^{8} such that their coordinates are the lengths of the eight edges of C(A) or, equivalently, the numbers in the chain code of C(A). More formally, for i = 1, 2, 3, 4, let n_{i}^{max}(A) and n_{i}^{min}(A) be the number of points x ∈ C(A) such that, respectively, x_{i} = x_{i}^{max}(A) and x_{i} = x_{i}^{min}(A), that is, n_{i}^{max}(A) = |{x ∈ C(A) : x_{i} = x_{i}^{max}(A)}| and n_{i}^{min}(A) = |{x ∈ C(A) : x_{i} = x_{i}^{min}(A)}|. The projection vector of A is the vector ν(A) ∈ ZZ^{8}, given by, for i = 1, 2, 3, 4, ν_{i}(A) = n_{i}^{max}(A) − 1 and ν_{i+4}(A) = n_{i}^{min}(A) − 1. For example,

the projection vector for A = [ \begin{matrix}1 & 1 & 0 \\ 0 & 1 & 1 \end{matrix} ] is ν(A) = (0, 1, 1, 0, 0, 1, 1, 0).

Note that the projection vector is independent of translation, that is, for any x ∈ ZZ^{2}, ν(A) = ν(A_{x}).

Proposition 8 Let A, X, Y ⊂ ZZ^{2}. If A = X ⊕ Y, then ν(A) = ν(X) + ν(Y).

The ideas for proving this proposition are given in [4].

Corollary 9 If A = B_{1} ⊕ B_{2} ⊕ ⋯ ⊕ B_{k}, then ν(A) = ν(B_{1}) + ν(B_{2}) + ⋯ + ν(B_{k}).

Corollary 10 Let A ⊂ ZZ^{2} and B_{1}, B_{2}, ⋯, B_{k} be a sequence of subsets of the elementary square, if ν(A) = ν(B_{1}) + ν(B_{2}) + ⋯ + ν(B_{k}) and A ≠ B_{1} ⊕ B_{2} ⊕ ⋯ ⊕ B_{k}, then B_{1}, B_{2}, ⋯, B_{k} is an unfeasible sequence for the decomposition of A.

Given two subsets X and Y, we denote ν(X) ≤ ν(Y) iff, for any i ∈ {1, …, 8}, ν_{i}(X) ≤ ν_{i}(Y).

Corollary 11 Let A ⊂ ZZ^{2} and B_{1}, B_{2}, ⋯, B_{k} be a sequence of subsets of the elementary square, if ν(A) ≱ ν(B_{1}) + ν(B_{2}) + ⋯ + ν(B_{k}), then B_{1}, B_{2}, ⋯, B_{k} is an unfeasible sequence for the decomposition of A.

Corollary 12 Let A ⊂ ZZ^{2} and B_{1}, B_{2}, ⋯, B_{m} be the invariant sequence of A. If ν(A) ≰ ν(B_{1}) + ν(B_{2}) + ⋯ + ν(B_{m}), then A is not decomposable.

3. Combinatorial Optimization techniques

The aim of optimization problems is to maximize or minimize a given function over a certain domain. Combinatorial Optimization problems are characterized by the case in which this domain is finite and its elements can be “easily” generated. The difficulty to solve such problems usually stays in the cardinality of the domain. It is not possible to enumerate, in a reasonable time, all possibilities in order to find the best solution.

A strategy to partially enumerate the solutions of a combinatorial optimization problem is known as the branch and bound method [7, p. 40]. The idea is to partition the set of possible solutions into some small sets, originating independent subproblems. When one solves all these subproblems, the best solution found is the optimal solution of the original problem. This partition can be viewed as the construction of a rooted tree, whose nodes correspond to partial solutions, and the root node to the original problem. The search in the tree can be improved if some lower and upper bounds on the value of the optimal solution are known. For example, if it is known a feasible solution whose value is say 25, and we are in a node with lower bound 28, this branch can be pruned, since the best possible solution in this branch cannot be better than a known solution. The art of the strategy is to find good bounds in order to avoid visiting all possibilities.

There are many different ways to visit the nodes of a tree. In the algorithm presented in Section 4, we use the depth first search [7, p. 39], that is described in the following. In the beginning all nodes are marked as “unvisited”. In the first iteration, some arbitrary node is selected (if the tree has a root, this is the selected node) and marked as “visited”. Its neighbors are marked “reachable”, and pushed in a stack. In an arbitrary iteration, if the stack is empty, the search stops. Otherwise, the top element is pop up and marked as “visited”. The still “unvisited” neighbors of this node are marked as “reachable”, pushed in the stack, and a new iteration begins.

4. Search of Optimal Decompositions

In this section, we describe a combinatorial optimization procedure for the search of decompositions of arbitrarily shaped SE’s by subsets of the elementary square. If there exist decompositions for a SE, it gives an optimal sequence of subsets, otherwise, it proves that the SE is not decomposable.

In order to test if a shape sequence represents a given SE, we use a classical matching algorithm. When there is a matching between the SE and the subset generated from the shape sequence, the translation sequence is immediately determined.

The hard part of the problem is the determination of the shape sequence. The formulation adopted is the branch and bound search, described in the last section, in a tree that represents the space of all possible sequences of subsets. In the search we look for feasible sequences, pruning unfeasible ones.

Corollary 2 and Proposition 4 show that the search space for determining decompositions of a SE A can be restricted to the sequence of invariants of A. Thus, given a SE A, we will define a tree that represents the space of all possible subsequences of the invariant sequence of A.

Let B_{1}, B_{2}, …, B_{m} be the invariant sequence of a SE A. The decomposition tree of A is a tree such that: i - the root is the unitary set that contains the origin; ii - a generic node is either the root or a subset Y = B_{i_{1}} ⊕ B_{i_{2}} ⊕ ⋯ ⊕ B_{i_{j}}, where B_{i_{1}}, B_{i_{2}}, …, B_{i_{j}} is a nomempty subsequence of the invariant sequence of A. iii - the direct descendents of a node Y = B_{i_{1}} ⊕ B_{i_{2}} ⊕ ⋯ ⊕ B_{i_{j}} are Y ⊕ B_{i_{j}+1}, Y ⊕ B_{i_{j}+2}, …, Y ⊕ B_{m}. iv - the edge that joins the node Y = B_{i_{1}} ⊕ B_{i_{2}} ⊕ ⋯ ⊕ B_{i_{j}} and its direct descendent Y ⊕ B_{k} is labeled B_{k} (k = i_{j} + 1, i_{j} + 2, …, m).

The level of a node X in a decomposition tree is the length of the path from the root to X and it is denoted by level(X).

Note that the labels of the edges on the path from the root to a node X forms a subsequence of the invariant sequence and that the subset X is given by Minkowski addition of the elements of this subsequence. So, if the node X is equal to A, then A has a decomposition given by the sequence defined by the path from the root to this node.

Moreover, observe that the nodes of the decomposition tree of A correspond to all possible 2^{m} subsequences of the invariant sequence of A. Our objective is to find the node A with minimum level. If there is no node A, then the SE is not decomposable.

A crucial task of the procedure is an adequate pruning of the decomposition tree. Corollary 12 gives a necessary condition for the existence of a decomposition for A.

Thus, we should begin the search just if this condition is satisfied. Proposition 6 fixes the lower bound for a decomposition of A, while Corollaries 10 and 11 guaranty important pruning, since they detect unfeasible sequences.

Under this context, after setting the lower bound (with the value given by Proposition 6) and the initial upper bound (with infinite), the search dynamics goes on. In the first iteration, the root is selected to be visited. In an arbitrary iteration, when the node X is being visited: i - if X satisfies the constraints of Corollary 11 (i.e., ν(A) ≱ ν(X)), or if X is not an invariant of A, or if level(X) is greater or equal to the current upper bound, then this node is pruned; ii - if level(X) is greater or equal to the lower bound and the constraint of Corollary 9 (i.e., ν(X) = ν(A)) is satisfied: ii.a - if X = A, then the path from the root to X is the present optimal solution and the upper bound is updated to level(X), additionally, if level(X) is equal to the lower bound, then X is the optimal solution; ii.b - if X ≠ A, by Corollary 10, the node X can be pruned; iii - if level(X) is less than the lower bound or the constraint of Corollary 9 is not satisfied (i.e., ν(X) ≠ ν(A)), then the node X is an invariant of A, ν(X) ≤ ν(A) and level(X) is less than the current upper bound, so the next iteration visits one of its unvisited descendents.

Depending on the order chosen to construct the invariant sequence, different heuristic search procedures arise. We have sorted the invariant sequence in decreasing order, according to the sum of the coordinates of the orthogonal projections of each subset.

5. Experimental Results

We have applied our search procedure for the decomposition of digital disks, convex sets, decomposable simply connected sets and decomposable sets that have holes.

The digital disk of radius r, centered at the origin, is the subset of ZZ^{2} given by D(r) = {(x, y) ∈ ZZ^{2} : x^{2} + y^{2} ≤ r^{2}}. For the disks of radius 2 and 4, the procedure spent less than 1 second to find the first decomposition, that were optimal. The disk D(3) and all disks of radius 5 to 16 are not decomposable. The time spent for detecting the non-decomposability of disks D(3), D(5) and D(6) was less than 1 second and for D(7) to D(16) was increasing from 2 seconds to 4 hours.

For convex sets the first solution found was always the optimal one. The time spent for finding the optimal decomposition of subsets of the 7 × 7 to 61 × 61 square were increasing from 0.6 to 50 seconds.

For decomposable simply connect sets the time spent, for example, for the first solution of a subset of the 46 × 46 square was 12 hours. Observing some experiments of this kind, we note that the distance between the optimal and the first solution was at most 1.

For decomposable sets that have holes the time spent for SE’s that are subsets of the 9 × 9, 17 × 17, 19 × 19, 27 × 27, 41 × 41 and 61 × 61 square were less than 22 seconds. All the first solutions were optimal.

6. Conclusion

We have presented a combinatorial optimization approach for finding optimal sequential decompositions of erosions and dilations by arbitrarily shaped SE’s.The procedure proposed gives an optimal decomposition or proves that the operator considered is not decomposable.

The central idea of the method is walking intelligently on the space of possible decompositions, supported by algebraic and geometrical properties of the SE’s.

We have presented a description of the search procedure as well as the Mathematical properties that justify the pruning. As should be expected, the complexity of the search procedure is greatly reduced for convex SE’s.

Several examples illustrate the capacity of the approach for finding decompositions of complex SE’s in reasonable processing times (compared with the size of the search space). These examples also show that the distance between the optimal and the first solution obtained usually is small. Therefore, in practical applications, it may be a good heuristic to stop when the first solution is obtained.

The model proposed is flexible and open for the addition of other pruning strategies, that should be a consequence of a deeper theoretical understanding of the problem.

The symbiosis between Mathematical Morphology and Combinatorial Optimization seems a powerful approach to solve hard problems as the optimization of the representation of morphological operators in the context of the automatic programming of morphological machines. The problem studied here is just a point in this universe.

References

1. J. Barrera, E. R. Dougherty, and N. S. Tomita. Automatic programming of binary morphological machines by design of statistically optimal operators in the context of computational learning theory. Journal of Electronic Imaging, 6(1):54–67, January 1997.

2. J. Barrera and R. F. Hashimoto. Compact representation of w-operators. Submited to International Symposium on Mathematical Morphology’ 98.

3. J. Barrera and G. P. Salas. Set Operations on Collections of Closed Intervals and their Applications to the Automatic Programming of Morphological Machines. Journal of Eletronic Imaging, 5(3):335–352, July 1996.

4. P. K. Ghosh. Vision, geometry, and minkowski operators. Contemporary Mathematics, 119:63– 83, 1991.

5. P. A. Maragos. A Unified Theory of Translation-invariant Systems with Applications to Morphological Analysis and Coding of Images. PhD thesis, School of Elect. Eng. - Georgia Inst. Tech., 1985.

6. H. Park and R. T. Chin. Decomposition of Arbitrarily Shaped Morphological Structuring Elements. IEEE Transactions on Pattern Analysis and Machine Intelligence, 17(1):2–15, January 1995.

7. E. Rich and K. Knight. Artificial Intelligence. McGraw-Hill, USA, second edition, 1991.

8. C. H. Richardson and R. Schafer. A lower bound for structuring element decompositions. IEEE Transactions on Pattern Analysis and Machine Intelligence, 13(4):365–369, April 1991.

9. J. Serra. Image Analysis and Mathematical Morphology. Academic Press, 1982.

10. J. Xu. Decomposition of Convex Polygonal Morphological Structuring Elements into Neighborhood Subsets. IEEE Transactions on Pattern Analysis and Machine Intelligence, 13(2):153– 162, February 1991.