Search code examples
pythonconditional-statementspolygongeopandas

Condition to harvest gpd.column values


How do I create a condition that will only harvest ['grnd', 'ot', 'rf'] values when the columns exist? i.e.: harvest data from these columns if present else continue.

import geopandas as gpd
import matplotlib.pyplot as plt

polys = gpd.GeoSeries([Polygon([(0,0), (2,0), (2, 1.5), (2,2), (0,2)]),
                     Polygon([(0,2), (2,2), (2,4), (0,4)]),
                     Polygon([(2,0), (5,0), (5,1.5), (2,1.5)]),
                     Polygon([(3,3), (5,3), (5,5), (3,5)]),
                       Polygon([(6,1), (8, 1), (8, 3), (6, 3)], 
                                     [[(6.5, 1.5), (7.5, 1.5), (7.5, 2.5), (6.5, 2.5)][::-1]]
                                            )])

fp = gpd.GeoDataFrame({'geometry': polys, 'name': ['a', 'b', 'c', 'd', 'e'],
                      'grnd': [25, 25, 25, 25, 25],
                      'rf': [29, 35, 26, 31, 28],
                      'ot': [np.nan, 31, 29, 32, 30]})

fig, ax = plt.subplots(figsize=(5, 5))
fp.plot(ax=ax, alpha=0.3, cmap='tab10', edgecolor='k',)
fp.apply(lambda x: ax.annotate(text=x['name'], xy=x.geometry.centroid.coords[0], ha='center'), axis=1)
plt.show()

enter image description here

We highlight polygons that share an edge (adjacent / touching) and the attributes of those polygons.

def extract_exterior_sides(polygon):
    if polygon.interiors:
        return list(map(LineString, zip(polygon.exterior.coords[:-1], polygon.exterior.coords[1:])))
    else:
        return list(map(LineString, zip(polygon.boundary.coords[:-1], polygon.boundary.coords[1:])))

# Extract all sides of all polygons (excluding sides of holes)
sides = fp.geometry.apply(lambda x: extract_exterior_sides(x)).explode()

with warnings.catch_warnings():
    warnings.simplefilter('ignore')
    result = {
        poly_name: {
            'sides': [
                side if poly_sides.index[i] % 2 == 0 else side.reverse() 
                for i, side in enumerate(poly_sides)
            ],
            'edges': [
                sorted([
                    i 
                    for i in fp.loc[
                        (fp.geometry.touches(line)) & (fp.geometry.intersection(line).length > 0), 
                        ['grnd', 'ot', 'rf']
                    ].values.tolist()[0] 
                    if not np.isnan(i)
                ])
                if len(fp.loc[(fp.geometry.touches(line)) & (fp.geometry.intersection(line).length > 0)]) == 1
                else sorted([
                    i 
                    for j in fp.loc[
                        (fp.geometry.touches(line)) & (fp.geometry.intersection(line).length > 0), 
                        ['grnd', 'ot', 'rf']
                    ].values.tolist() 
                    for i in j 
                    if not np.isnan(i)
                ])
                for line in [
                    side if poly_sides.index[i] % 2 == 0 else side.reverse() 
                    for i, side in enumerate(poly_sides)
                ]
            ]
        }
        for poly_name, poly_sides in sides.groupby(fp.name)
    }

print(result)

which yields {'a': {'sides': [<LINESTRING (0 0, 2 0)>, <LINESTRING (2 0, 2 1.5)>, <LINESTRING (2 1.5, 2 2)>, <LINESTRING (2 2, 0 2)>, <LINESTRING (0 2, 0 0)>], 'edges': [[25.0, 29.0], [25.0, 25.0, 26.0, 29.0, 29.0], [25.0, 29.0], [25.0, 25.0, 29.0, 31.0, 35.0], [25.0, 29.0]]}, 'b': {'sides': [<LINESTRING (2 2, 0 2)>, <LINESTRING (2 4, 2 2)>, <LINESTRING (0 4, 2 4)>, <LINESTRING (0 2, 0 4)>], 'edges': [[25.0, 25.0, 29.0, 31.0, 35.0], [25.0, 31.0, 35.0], [25.0, 31.0, 35.0], [25.0, 31.0, 35.0]]}, 'c': {'sides': [<LINESTRING (2 0, 5 0)>, <LINESTRING (5 0, 5 1.5)>, <LINESTRING (5 1.5, 2 1.5)>, <LINESTRING (2 1.5, 2 0)>], 'edges': [[25.0, 26.0, 29.0], [25.0, 26.0, 29.0], [25.0, 26.0, 29.0], [25.0, 25.0, 26.0, 29.0, 29.0]]}, 'd': {'sides': [<LINESTRING (5 3, 3 3)>, <LINESTRING (5 5, 5 3)>, <LINESTRING (3 5, 5 5)>, <LINESTRING (3 3, 3 5)>], 'edges': [[25.0, 31.0, 32.0], [25.0, 31.0, 32.0], [25.0, 31.0, 32.0], [25.0, 31.0, 32.0]]}, 'e': {'sides': [<LINESTRING (6 1, 8 1)>, <LINESTRING (8 1, 8 3)>, <LINESTRING (8 3, 6 3)>, <LINESTRING (6 3, 6 1)>], 'edges': [[25.0, 28.0, 30.0], [25.0, 28.0, 30.0], [25.0, 28.0, 30.0], [25.0, 28.0, 30.0]]}}

So; if the fp GeoDataFrame had no 'grnd': [25, 25, 25, 25, 25] column; how would the code execute successfully?


Solution

  • The column names that you select from the data frame don't have to be hard coded -- they're just an ordinary Python list. So you can do this by building a list of columns to work with at run time. You could do this by adding

    cols = [c for c in ['grnd', 'ot', 'rf'] if c in fp.columns]
    

    before the result = line. This will give you the subset of interesting columns that are actually present in fp.

    Then replace the two mentions of ['grnd', 'ot', 'rf'] later in the code with cols.