# Creating mask of geospatial data with maximum side-length(triangulation)

Earlier I asked a question, about Python GeoPandas crops result or sets boundary to data

And I still cant understand how to work with geospatial data.

I generate coordinates and create a plot with this code:

``````import csv

import numpy
from matplotlib import pyplot
from matplotlib.tri import Triangulation

result = list()

with open(path) as csv_handler:

for item in line[0].split(';'):
item = float(item)

if item == -999.0:
result.append(numpy.nan)
elif item == -999.7999877929688:
result.append(numpy.nan)
else:
result.append(item)

return result

def nan_equal(actual):
try:
numpy.testing.assert_equal(actual, numpy.nan)
except AssertionError:
return False

return True

def clear(longitudes, latitudes, values):
if len(longitudes) == len(latitudes) == len(values):
result_longitudes = list()
result_latitudes = list()
result_values = list()

for i in range(len(longitudes)):
if not nan_equal(longitudes[i]):
if not nan_equal(latitudes[i]):
if not nan_equal(values[i]):
result_longitudes.append(longitudes[i])
result_latitudes.append(latitudes[i])
result_values.append(values[i])

return result_longitudes, result_latitudes, result_values
else:
raise ValueError

def main():
size = 5000

longitudes = numpy.random.rand(size)
latitudes = numpy.random.rand(size)
values = numpy.random.uniform(-1000.0, 1000.0, size)

triangulation1 = Triangulation(longitudes, latitudes)
triangulation2 = Triangulation(longitudes, latitudes)
triangulation3 = Triangulation(longitudes, latitudes)

triangles = triangulation.triangles

xtri = longitudes[triangles] - numpy.roll(longitudes[triangles], 1, axis=1)
ytri = latitudes[triangles] - numpy.roll(latitudes[triangles], 1, axis=1)
maxi = numpy.max(numpy.sqrt(xtri ** 2 + ytri ** 2), axis=1)

figure, (axes1, axes2, axes3) = pyplot.subplots(ncols=3, figsize=(60, 20))

axes1.tricontourf(triangulation1, values)
axes1.scatter(longitudes, latitudes, s=5, color='black')

axes2.tricontourf(triangulation2, values)
axes2.scatter(longitudes, latitudes, s=5, color='black')

axes3.tricontourf(triangulation3, values)
axes3.scatter(longitudes, latitudes, s=5, color='black')

pyplot.show()

if __name__ == '__main__':
main()
``````

But if i change random coordinates generation to:

``````longitudes = numpy.random.uniform(-180.0, 180.0, size)
latitudes = numpy.random.uniform(-90.0, 90.0, size)
``````

Code doesnt work and write error message

``````ValueError: z array must not contain non-finite values within the triangulation
``````

Update

I found some code to create polygon by points, but I don't know how to convert it to triangular mask

``````def alpha_shape(points, alpha):
if len(points) < 4:
return geometry.MultiPoint(list(points)).convex_hull

def add_edge(edges, edge_points, coords, i, j):
if (i, j) in edges or (j, i) in edges:
return
edge_points.append(coords[[i, j]])

coords = numpy.array([point.coords[0] for point in points])

tri = Delaunay(coords)
edges = set()
edge_points = []
for ia, ib, ic in tri.vertices:
pa = coords[ia]
pb = coords[ib]
pc = coords[ic]

a = math.sqrt((pa[0] - pb[0]) ** 2 + (pa[1] - pb[1]) ** 2)
b = math.sqrt((pb[0] - pc[0]) ** 2 + (pb[1] - pc[1]) ** 2)
c = math.sqrt((pc[0] - pa[0]) ** 2 + (pc[1] - pa[1]) ** 2)

s = (a + b + c) / 2.0

area = math.sqrt(s * (s - a) * (s - b) * (s - c))
circum_r = a * b * c / (4.0 * area)

if circum_r < 1.0 / alpha:

m = geometry.MultiLineString(edge_points)
triangles = list(polygonize(m))