I have a polygon GeoDataFrame 'Polygons' with attribute
"AREA" (total polygon area in ha) and another GeoDataFrame 'Small_polygon' with smaller polygons within the 'Polygons'.
How do I calculate smaller polygon area within large polygons?
What I have already is Python code where I loop through each polygon in 'Polygons' and find corresponding smaller polygons within the 'Polygons' using the index.
After that I need to find each smaller polygon area within a large polygon.
For example, here is a polygon with
Polygons['AREA'] = 100 ha and the result is calculated area (ha) for each small polygon.
I found some solutions in ArcGIS and QGIS, however I would like this to be done using Python.
I do not need to store the result anywhere, I will use them for further calculations within a loop, so this is kind of middle stage of the analyses. The further steps would be using the 'Some_points' layer calculate each smaller polygon average and max diameter, height and so on.
import pandas as pd import geopandas as gpd Polygons = gpd.read_file(r'E:\...\Polygons.shp') Small_polygon = gpd.read_file(r'E:\...\Small_polygon.shp') some_points = gpd.read_file(r'E:\...\points.shp') print(Polygons.crs == Small_polygon.crs == some_points.crs) True # define the ID for each polygon Polygons['ID'] = Polygons.index id = Polygons['ID'] Polygons.drop(labels=['ID'], axis=1, inplace=True) Polygons.insert(0, 'ID', id) small_poly_join = gpd.sjoin(Small_polygon, Polygons)[['ID', 'small_id', 'geometry']] points_join = gpd.sjoin(small_poly_join, some_points)[['ID', 'small_id', 'height', 'diameter']] for i in Polygons.index: df_1 = small_poly_join[small_poly_join['FID'] == i] # this gives a wrong result area = (df_1['geometry'].area/ 10**6) # km2 area = (df_1['geometry'].area/10000 # ha for n in df_1.index: df_2 = points_join[points_join['id'] == n] max_height = df_2['height'].max() ave_height = df_2['height'].mean() max_diam = df_2['diameter'].max() ave_diam = df_2['diameter'].mean()
I had some problems with the area calculation that I found in other question, for instance,
area = df_1['geometry'].area / 10000 gives me a wrong result compared to total polygon area. I also checked the CRS and projection which are equal to all used shapefiles.