I am trying to find the spatial relationship between two polygon geoDFs as follows:

  1. Find geoDF1.geom within geoDF2.geom else (not within)
  2. Find geoDF1.geom with largest intersection fromgeoDF2.geom else (no intersection)
  3. Find geoDF1.geom with smallest distance from geoDF2.geom
  4. Update geoDF1 with selected column from from geoDF2

I wrote up a function which seems to works fine, but I am just not sure this is the most effective way to do it. I am especially worried about step #3 which right now checks the distance from each one of geoDF2 rows, which can become memory intensive if it is a large table. Perhaps I should create a buffer that limits the distance checked?

import geopandas as gpd

### FUNCTION that checks spatial relationship between two geodata frames - df1, df2:
# Order of operations: Within -> largest intersection -> smallest distance
# Checks relationship of df1.geom (left) to df2.geom
# Updates df1 with new col from df2:: df2_col
# df1_id = unique id

def relate(df1, df1_id, df2, df2_col):

    # join df1_col for all df1_geom within df2_geom
    join_df = gpd.sjoin(df1, df2, how='left', op="within")
    ##sort by col_1
    join_df = join_df.sort_values(df1_id)

    # get rows that are not within
    remaining_df = join_df.loc[join_df[df2_col].isna()].reset_index(drop=True).drop(columns='index_right')

    ##checks to see relations of df1 to df2
    intersects = pd.DataFrame(columns=['id', 'type'])
    # check for each geom1 not within geom2
    for df1_geom, id in zip(remaining_df.geom, remaining_df[df1_id]):
        size = 0
        for df2_geom, type in zip(df2.geom, df2[df2_col]):
            # see size of intersection , if larger (or first) than previous intersection, set 'sel_type' to df2_col to associate with df2 (largest intersection)
                new_size = df1_geom.intersection(df2_geom).area
                if new_size > size:
                    size = new_size
                    sel_type = type
                    # print("intersects" + str(id) + sug)
                print("couldn't intersect for df1 with id: " + str(id))
        # if intersection exists add to lists of df1_id, df2_col
        if size > 0:
            intersects = intersects.append({'id': id, 'type': sel_type}, ignore_index=True)

    # sort by df1_id
    ints = intersects.sort_values('id').set_index('id', drop=True)

    # set type for intersecting geoms
    join_df.loc[join_df[df1_id].isin(ints.index), df2_col] = list(ints['type'])

    # #get df1 geom, who are not within and do not intersect with df2_geom
    remaining_df2 = join_df.loc[join_df[df2_col].isna()].reset_index(drop=True).drop(columns='index_right')

    ##!!## JOIN based on distance
    distance = pd.DataFrame(columns=['id', 'type'])
    # check for each df1 with no df2_col
    for geom, id in zip(remaining_df2.geom, remaining_df2[df1_id]):
        size = None
        for df2_geom, type in zip(df2.geom, df2[df2_col]):
            # see size of distance from df2_geom, if smaller (or first) than previous distance, set 'sel_type' to associate with df1_geom (minimal distance)
                insize = geom.distance(df2_geom)
                if size == None or insize < size:
                    size = insize
                    sel_type = type
                print("couldn't calculate distance for df1 with id: " + str(id))
        if size != None:
            distance = distance.append({'id': id, 'type': sel_type}, ignore_index=True)

        # sort by df1 id
    dis = distance.sort_values('id').set_index('id', drop=True)

    # set df2_col type for remaining df1 geoms by minimal distance
    join_df.loc[join_df[df1_id].isin(dis.index), df2_col] = list(dis['type'])

    ##UPDATE columns on original df1 with col to join from df2
    df1[df2_col] = [join_df.loc[join_df[df1_id==x], df2_col] for x in df1[df1_id]]

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.