# Create a structuringselement in the form of a line with a certain degree and length

Within Matlab there exist a function to create a structuringselement in the form of a line.
SE = strel('line',len,deg) creates a linear structuring element that is symmetric with respect to the neighborhood center, with approximate length len and angle deg. For a couple of morphological operations I need to filter out horizontal-, vertical and other angles of lines. Within Python's CV2 a similiar function does not exist. Therefore, I wanted to recreate this function in Python

My code:

import math
import numpy as np

# bresenham function is the accepted answer of SO's post https://stackoverflow.com/questions/23930274/list-of-coordinates-between-irregular-points-in-python
def bresenham(x0, y0, x1, y1):
points = []
dx = abs(x1 - x0)
dy = abs(y1 - y0)
x, y = x0, y0
sx = -1 if x0 > x1 else 1
sy = -1 if y0 > y1 else 1
if dx > dy:
err = dx / 2.0
while x != x1:
points.append((x, y))
err -= dy
if err < 0:
y += sy
err += dx
x += sx
else:
err = dy / 2.0
while y != y1:
points.append((x, y))
err -= dx
if err < 0:
x += sx
err += dy
y += sy
points.append((x, y))

return points

def strel_line(length, degrees):
if length >= 1:
theta = degrees * np.pi / 180
x = round((length - 1) / 2 * np.cos(theta))
y = -round((length - 1) / 2 * np.sin(theta))
points = bresenham(-x, -y, x, y)
points_x = [point[0] for point in points]
points_y = [point[1] for point in points]
n_rows = int(2 * max([abs(point_y) for point_y in points_y]) + 1)
n_columns = int(2 * max([abs(point_x) for point_x in points_x]) + 1)
strel = np.zeros((n_rows, n_columns))
rows = ([point_y + max([abs(point_y) for point_y in points_y]) for point_y in points_y])
columns = ([point_x + max([abs(point_x) for point_x in points_x]) for point_x in points_x])
idx = []
for x in zip(rows, columns):
idx.append(np.ravel_multi_index((int(x[0]), int(x[1])), (n_rows, n_columns)))
strel.reshape(-1)[idx] = 1

return strel

if __name__=='__main__':
strel = strel_line(15, 105)
print(strel)
print(strel.shape)


Output:

[[1. 0. 0. 0. 0.]
[1. 0. 0. 0. 0.]
[0. 1. 0. 0. 0.]
[0. 1. 0. 0. 0.]
[0. 1. 0. 0. 0.]
[0. 1. 0. 0. 0.]
[0. 0. 1. 0. 0.]
[0. 0. 1. 0. 0.]
[0. 0. 1. 0. 0.]
[0. 0. 0. 1. 0.]
[0. 0. 0. 1. 0.]
[0. 0. 0. 1. 0.]
[0. 0. 0. 1. 0.]
[0. 0. 0. 0. 1.]
[0. 0. 0. 0. 1.]]


Do you see any ways to write my function more efficiently?

## Type hints

def bresenham(x0, y0, x1, y1):


should be, at a guess,

def bresenham(x0: float, y0: float, x1: float, y1: float) -> List[float]:


or maybe int, but you get the idea.

## Type promotion

You don't need decimals here:

err = dx / 2.0


/2 will accomplish the same thing.

## Termination

What happens if you call bresenham(0, 0, 9.5, 9.5)? I suspect that it will iterate infinitely and hang, due to

  while x != x1:


One way to fix this, if your algorithm permits, is - rather than making an sx dependent on the comparison of x0 and x1, simply swap x0 and x1 if they're in the wrong order. Then, rather than !=, use <, and you don't have to calculate sx.

## Use a generator

For some use cases this will dramatically decrease memory consumption. Instead of building up points in memory, just yield whenever you have a point. Allow the caller to decide whether they want to cast to a list, a tuple, or just iterate over the points.

## Int casting

  n_rows = int(2 * max([abs(point_y) for point_y in points_y]) + 1)


is overambitious in application of the int cast. The multiplication and addition don't need to be done in floating-point math, so you can instead do

  n_rows = 2 * int(max(abs(point_y) for point_y in points_y)) + 1


Also note that you shouldn't pass a list to max - just pass the generator directly.