using shapely str tree on coordinates to find points within x meters

shapely_rtree.py

from shapely.geometry import MultiPoint, Point                                                        
from shapely.ops import transform                                                         
from shapely.strtree import STRtree                                                                                                                                                                                
                                                                                          
import pandas as pd                                                                       
import pyproj                                                                                                                                                                                                      
                                                                                                                                                                                                                   
                                                                                                                                                                                                                   
def main(point, coordinates):                                                                                                                                                                                      
    points = [Point(longitude, latitude) for latitude, longitude in coordinates]                                                                                                                                   
    mp = MultiPoint(points)                                                                                                                                                                                        
    myProj = pyproj.Proj(                                                                                                                                                                                          
        "+proj=utm +zone=23K, +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs"                                                                                                                                  
    )                                                                                                                                                                                                              
    transformed_points = transform(myProj, mp)                                                                                                                                                                     
    transformed_point = transform(myProj, point)                                                                                                                                                                   
    max_distance_in_meters = 804.672  # about 1/2 a mile                                                                                                                                                           
    tree = STRtree(transformed_points)                                                                                                                                                                             
    close_points = tree.query(transformed_point.buffer(max_distance_in_meters))                                                                                                                                    
    return [i for i, p in enumerate(transformed_points) if p in close_points]                                                                                                                                      
                                                                                                                                                                                                                   
                                                                                                                                                                                                                   
if __name__ == "__main__":                                                                                                                                                                                         
    lng = -74.0173461                                                                                                                                                                                              
    lat = 40.7044542                                                                                                                                                                                               
    point = Point(lng, lat)                                                                                                                                                                                        
                                                                                                                                                                                                                   
    df = pd.read_csv("./points.csv")                                                                                                                                                                               
    coordinates = df[["latitude", "longitude"]].values                                                                                                                                                             
    indices = main(point, coordinates)                                                                                                                                                                             
    print(df.iloc[indices])