import pandas as pd
import matplotlib.pyplot as plt
import geopandas as gpd
import folium
import contextily as cx
import rtree
from zlib import crc32
import hashlib
from shapely.geometry import Point, LineString, Polygon
import numpy as np
from scipy.spatial import cKDTree
from shapely.geometry import Point
from haversine import Unit
from geopy.distance import distance

## Importing our DataFrames

gisfilepath = "/Users/jnapolitano/Projects/data/energy/non-active-wells.geojson"


wells_df = gpd.read_file(gisfilepath)

wells_df = wells_df.to_crs(epsg=3857)

## Importing our DataFrames

gisfilepath = "/Users/jnapolitano/Projects/data/energy/Liquified_Natural_Gas_Import_Exports_and_Terminals.geojson"


terminal_df = gpd.read_file(gisfilepath)

terminal_df = terminal_df.to_crs(epsg=3857)

terminal_df.drop(terminal_df[terminal_df['STATUS'] == 'SUSPENDED'].index, inplace = True)
terminal_df.rename(columns={"NAME": "TERMINAL_NAME"})
terminal_df['TERMINAL_GEO'] = terminal_df['geometry'].copy()
terminal_df.columns

terminal_map =terminal_df.explore(
    column="TYPE", # make choropleth based on "PORT_NAME" column
     popup=True, # show all values in popup (on click)
     tiles="Stamen Terrain", # use "CartoDB positron" tiles
     cmap='Set1', # use "Set1" matplotlib colormap
     #style_kwds=dict(color="black"),
     marker_kwds= dict(radius=6),
     #tooltip=['NAICS_DESC','REGION', 'COMMODITY' ],
     legend =True, # use black outline)
     categorical=True,
    )


terminal_map

def ckdnearest(gdA, gdB):

    nA = np.array(list(gdA.geometry.apply(lambda x: (x.x, x.y))))
    nB = np.array(list(gdB.geometry.apply(lambda x: (x.x, x.y))))
    btree = cKDTree(nB)
    dist, idx = btree.query(nA, k=1)
    gdB_nearest = gdB.iloc[idx].drop(columns="geometry").reset_index(drop=True)
    gdf = pd.concat(
        [
            gdA.reset_index(drop=True),
            gdB_nearest,
            pd.Series(dist, name='dist')
        ], 
        axis=1)

    return gdf

c = ckdnearest(wells_df, terminal_df)

c.describe()

nearest_wells_df= wells_df.sjoin_nearest(terminal_df, distance_col="distance_euclidian")

nearest_wells_df.describe()

#df.geopy.distance.distance(coords_1, coords_2).km
#df.apply(lambda row: distance(row['point'], row['point_next']).km if row['point_next'] is not None else float('nan'), axis=1)
# Thanks to https://stackoverflow.com/questions/55909305/using-geopy-in-a-dataframe-to-get-distances

nearest_wells_df['true_distance_km'] = nearest_wells_df.apply(lambda row: distance((row['LATITUDE_left'], row['LONGITUDE_left']), (row['LATITUDE_right'], row['LONGITUDE_right'])).km if row['geometry'] is not None else float('nan'), axis=1)

nearest_wells_df.describe()

filtered_wells = nearest_wells_df.loc[nearest_wells_df['true_distance_km'] < 50].copy()

filtered_wells.describe()

filtered_wells.explore(
    column="STATUS_left", # make choropleth based on "PORT_NAME" column
     popup=True, # show all values in popup (on click)
     tiles="Stamen Terrain", # use "CartoDB positron" tiles
     cmap='Set1', # use "Set1" matplotlib colormap
     #style_kwds=dict(color="black"),
     marker_kwds= dict(radius=6),
     #tooltip=['NAICS_DESC','REGION', 'COMMODITY' ],
     legend =True, # use black outline)
     categorical=True,)