0

I have defined two space dimesions ( x and z ) and I was able to manually "draw" an object to use it as a boolen for solving an equation. I defined it as it follows:

A = np.zeros((nz,nx))
object = np.ones_like(A)
object[ int(5/dz):int(10/dz) , int(5/dx):int(10/dz) ] = 2
object = object == 2

By doing that I can define an square 5x10 in z dimesion and 5x10 in x dimesion , and apply the algorythim which understands this as an area , I think. But when it comes to draw complex areas it ends up being hard doing it by little squares and rectangles.

So I want to automatize an area generation by mouse clicking and I want to be able to use this area as a boolean.

I was able to draw a polygon using:

import matplotlib.pyplot as plt
import numpy as np
from matplotlib.patches import Polygon

fig, ax = plt.subplots()

object = np.array(plt.ginput(n=-100,mouse_stop=2))
p = Polygon(object, alpha=0.5)
plt.gca().add_artist(p)


plt.draw()
plt.show()

But this outputs z and x coordinates of the vertices, and I tried to use it as boleean but I could'nt write it so that python uderstands it as the area defined by those points.

Is this problem easy to solve?

Tepe86
  • 1
  • 1
  • To clarify, what you want to do is to [rasterize](https://en.wikipedia.org/wiki/Rasterisation) an arbitrary polygon into a boolean mask? – Aleon Nov 14 '19 at 13:45
  • I want to be able to use the area of an arbitrary polygon as a name that I define (boolean), so I'm able to call diferente properties for this polygon (that should exist inside z and x dimesions). Like this example, for K being condutivity: `code`K = Kgeneral*np.ones_like(A) / K[object] = 100`code` 'object' being the area of my polygon, where I want to add a diferent property in regards to its surrouding. – Tepe86 Nov 14 '19 at 14:22
  • As I searched in google, I think it is exatly what you said, if I can rasterize the polygon, it would work because it would generate z and x coordinates as squares or rectangles that compose the polygon and I would be able to use this. – Tepe86 Nov 14 '19 at 17:56

1 Answers1

0

If you just want to calculate the area of a general polygon, you can use for example the Shapely python package like this:

import numpy as np
import matplotlib.pyplot as plt
from shapely.ops import Polygon
from matplotlib.patches import Polygon as PltPolygon

# Get the coordinate input
canvas_size = np.array([1, 1])
canvas_lim = np.array([[0, canvas_size[0]], [0, canvas_size[1]]])
fig, ax = plt.subplots()
plt.xlim(canvas_lim[0])
plt.ylim(canvas_lim[1])
ax.set_aspect("equal")
coordinates = np.array(plt.ginput(n=-100, mouse_stop=2))

# Use shapely.ops.Polygon to calculate the area
poly = Polygon(coordinates)

area = poly.area
print("The area is {} units^2".format(area))

# Draw the polygon
p = PltPolygon(coordinates, alpha=0.5)
ax.add_artist(p)
plt.show()

If you definitely need the mask, here's one way to rasterize it using numpy and matplotlib.path. For details see the comments in the code:

import numpy as np
import matplotlib.path as mpltPath
import matplotlib.pyplot as plt

# Define the limits of our polygon
canvas_desired_size = np.array([110, 100])
# The pixel size with which we calculate (number of points to consider)
# The higher this number, the more we have to calculate, but the
# closer the approximation will be
pixel_size = 0.1

# Cacluate the actual size of the canvas
num_pxiels = np.ceil(canvas_desired_size / pixel_size).astype(int)
canvas_actual_size = num_pxiels * pixel_size

# Let's create a grid where each pixel's value is it's position in our 2d image
x_coords = np.linspace(
    start=0,
    stop=canvas_actual_size[0],
    endpoint=False,
    num=canvas_desired_size[0] / pixel_size,
)
y_coords = np.linspace(
    start=0,
    stop=canvas_actual_size[1],
    endpoint=False,
    num=canvas_desired_size[1] / pixel_size,
)
# Since it makes more sense to check if the middle of the pixel is in the
# polygion, we shift everything with half pixel size
pixel_offset = pixel_size / 2
x_centers = x_coords + pixel_offset
y_centers = y_coords + pixel_offset

xx, yy = np.meshgrid(x_centers, y_centers, indexing="ij")

# Flatten our xx and yy matrixes to an N * 2 array, which contains
# every point in our grid
pixel_centers = np.array(
    list(zip(xx.flatten(), yy.flatten())), dtype=np.dtype("float64")
)

# Now prompt for the imput shape
canvas_lim = np.array([[0, canvas_actual_size[0]], [0, canvas_actual_size[1]]])
fig, ax = plt.subplots()
plt.xlim(canvas_lim[0])
plt.ylim(canvas_lim[1])
ax.set_aspect("equal")
shape_points = np.array(plt.ginput(n=-100, mouse_stop=2))

# Create a Path object
shape = mpltPath.Path(shape_points)
# Use Path.contains_points to calculate if each point is
# within our shape
shape_contains = shape.contains_points(pixel_centers)

# Reshape the result to be a matrix again
mask = np.reshape(shape_contains, num_pxiels)

# Calculate area
print(
    "The shape area is roughly {} units^2".format(
        np.sum(shape_contains) * pixel_size ** 2
    )
)

# Show the rasterized shape to confirm it looks correct
plt.imshow(np.transpose(mask), aspect="equal", origin="lower")
plt.xlim([0, num_pxiels[0]])
plt.ylim([0, num_pxiels[1]])
plt.show()

Alternatively, a simpler solution would be using your plot as an image and thresholding it to get a boolean mask. There should be plent of examples of how to do this on google.

Aleon
  • 311
  • 2
  • 10