Coverage for pyguymer3/geo/area.py: 66%
38 statements
« prev ^ index » next coverage.py v7.9.2, created at 2025-07-08 18:47 +0000
1#!/usr/bin/env python3
3# Define function ...
4def area(
5 shape,
6 /,
7 *,
8 eps = 1.0e-12,
9 level = 1,
10 nIter = 100,
11 onlyValid = False,
12 repair = False,
13):
14 """Find the area of a shape.
16 Parameters
17 ----------
18 shape : shapely.coords.CoordinateSequence, shapely.geometry.point.Point, shapely.geometry.multipoint.MultiPoint, shapely.geometry.polygon.LinearRing, shapely.geometry.linestring.LineString, shapely.geometry.multilinestring.MultiLineString, shapely.geometry.polygon.Polygon, shapely.geometry.multipolygon.MultiPolygon
19 the shape
20 eps : float, optional
21 the tolerance of the Vincenty formula iterations
22 level : int, optional
23 the number of levels to split the shape into
24 nIter : int, optional
25 the maximum number of iterations (particularly the Vincenty formula)
26 onlyValid : bool, optional
27 only add valid Polygons (checks for validity can take a while, if being
28 being called often)
29 repair : bool, optional
30 attempt to repair invalid Polygons
32 Returns
33 -------
34 area : float
35 the area (in metres-squared)
37 Notes
38 -----
39 Copyright 2017 Thomas Guymer [1]_
41 References
42 ----------
43 .. [1] PyGuymer3, https://github.com/Guymer/PyGuymer3
44 """
46 # Import special modules ...
47 try:
48 import shapely
49 import shapely.ops
50 except:
51 raise Exception("\"shapely\" is not installed; run \"pip install --user Shapely\"") from None
53 # Import sub-functions ...
54 from ._area import _area
55 from .extract_polys import extract_polys
57 # **************************************************************************
59 # Initialize total ...
60 tot = 0.0 # [m2]
62 # Loop over the Polygons in the shape ...
63 for shapePart in extract_polys(shape, onlyValid = onlyValid, repair = repair):
64 # Loop over the Polygons in the Voronoi diagram of the part of the shape ...
65 for voronoi in extract_polys(shapely.ops.voronoi_diagram(shapePart), onlyValid = onlyValid, repair = repair):
66 # Loop over the parts of the Polygon in the Voronoi diagram of the
67 # part of the shape which intersect the part of the shape ...
68 for voronoiPart in extract_polys(shapePart.intersection(voronoi), onlyValid = onlyValid, repair = repair):
69 # Loop over triangles within the Polygon ...
70 for triangle1 in extract_polys(shapely.ops.triangulate(voronoiPart), onlyValid = onlyValid, repair = repair):
71 # Check if the user wants this level of refinement ...
72 if level == 1:
73 # Increment the total and move on to the next one ...
74 tot += _area(triangle1, eps = eps, nIter = nIter) # [m2]
75 continue
77 # Loop over triangles within the triangle ...
78 for triangle2 in extract_polys(shapely.ops.triangulate(triangle1), onlyValid = onlyValid, repair = repair):
79 # Check if the user wants this level of refinement ...
80 if level == 2:
81 # Increment the total and move on to the next one ...
82 tot += _area(triangle2, eps = eps, nIter = nIter) # [m2]
83 continue
85 # Loop over triangles within the triangle ...
86 for triangle3 in extract_polys(shapely.ops.triangulate(triangle2), onlyValid = onlyValid, repair = repair):
87 # Check if the user wants this level of refinement ...
88 if level == 3:
89 # Increment the total and move on to the next
90 # one ...
91 tot += _area(triangle3, eps = eps, nIter = nIter) # [m2]
92 continue
94 # Loop over triangles within the triangle ...
95 for triangle4 in extract_polys(shapely.ops.triangulate(triangle3), onlyValid = onlyValid, repair = repair):
96 # Check if the user wants this level of
97 # refinement ...
98 if level == 4:
99 # Increment the total and move on to the
100 # next one ...
101 tot += _area(triangle4, eps = eps, nIter = nIter) # [m2]
102 continue
104 # Loop over triangles within the triangle ...
105 for triangle5 in extract_polys(shapely.ops.triangulate(triangle4), onlyValid = onlyValid, repair = repair):
106 # Check if the user wants this level of
107 # refinement ...
108 if level == 5:
109 # Increment the total and move on to the
110 # next one ...
111 tot += _area(triangle5, eps = eps, nIter = nIter) # [m2]
112 continue
114 # Loop over triangles within the triangle ...
115 for triangle6 in extract_polys(shapely.ops.triangulate(triangle5), onlyValid = onlyValid, repair = repair):
116 # Check if the user wants this level of
117 # refinement ...
118 if level == 6:
119 # Increment the total and move on to
120 # the next one ...
121 tot += _area(triangle6, eps = eps, nIter = nIter) # [m2]
122 continue
124 # Cry ...
125 raise Exception(f"\"level\" is too large (\"{level:,d}\")") from None
127 # Return total ...
128 return tot