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

Question

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


def read(path):
    result = list()

    with open(path) as csv_handler:
        csv_reader = csv.reader(csv_handler)

        for line in csv_reader:
            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)

    def set_mask(triangulation, alpha=0.4):
        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)

        triangulation.set_mask(maxi > alpha)

    set_mask(triangulation2, alpha=0.1)
    set_mask(triangulation3, alpha=0.3)

    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
        edges.add((i, j))
        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:
            add_edge(edges, edge_points, coords, ia, ib)
            add_edge(edges, edge_points, coords, ib, ic)
            add_edge(edges, edge_points, coords, ic, ia)

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

    concave_hull = cascaded_union(triangles)

    return concave_hull, edge_points