EPSG 4326 and EPSG 3857 of Web Mercator

buffer_ray.py

from pyproj import Transformer
from shapely.wkt import loads

transformer_4326_3857 = Transformer.from_crs("epsg:4326", "epsg:3857")
transformer_3857_4326 = Transformer.from_crs("epsg:3857", "epsg:4326")


def generate_buffer_ray(wkt_linestring, distance=15):
    lineString_4326 = loads(wkt_linestring)
    # Reproject LineString
    array_coords_3857 = []
    for point in list(lineString_4326.coords):
        x = point[0]
        y = point[1]
        x, y = transformer_4326_3857.transform(y, x)
        array_coords_3857.append((y, x))
    lineString_3857 = LineString(array_coords_3857)
    # # Buffer the lineString
    polygon_3857 = lineString_3857.buffer(distance)
    proj_buffer_points = []
    for point in list(polygon_3857.exterior.coords):
        x = point[0]
        y = point[1]
        x, y = transformer_3857_4326.transform(y, x)
        proj_buffer_points.append((y, x))
    polygon_4326 = Polygon(proj_buffer_points)
    bounds = polygon_4326.bounds
    bounds_poly = box(minx=bounds[0], miny=bounds[1], maxx=bounds[2], maxy=bounds[3])
    return bounds_poly.wkt

Difference between EPSG 4326 & 3857

EPSG: 4326 uses a coordinate system on the surface of a sphere or ellipsoid of reference.
WGS 84 - Earth as Geoid. -Mercator

EPSG: 3857 uses a coordinate system PROJECTED from the surface of the 
sphere. Earth as perfectly sphere. -Web Mercator

Think of it as this way:

EPSG 4326 uses a coordinate system the same as a GLOBE (curved surface). 
EPSG 3857 uses a coordinate system the same as a MAP (flat surface).


Convert GDAL to convert HARN NAD83 Washington North (EPSG of 2926) into Web mercator WGS84 (EPSG of 4326)

From NICK:
TileMill and Github are going to be using the same, Spherical Mercator. 
EPSG4326 is Plate Carree (Geographic), meaning that it is NOT yet transformed to a projection. 
You want your GeoJSON to be this. Leaflet likes this and does the projection to Web Mercator for you.

An interesting tid bit about 900913 is that it is GOOGLE in digits (like when you play with a calculator). 
They did this until EPSG made them a valid projection. 

web mercator NOT == spherical mercator 
google mercator == Web Mercator (Auxiliary Sphere)

Loading data to Mapbox Studio/TileMill has to be in WGS 84, Mercator EPSG:4326
Once those data gets uploaded, it converts into EPSG:3857 WebMercator.

distance.py

from shapely.geometry import LineString
from shapely.ops import transform
from functools import partial
import pyproj


def distance(point1, point2):
    line = LineString([(point1[0], point1[1]), (point2[0], point2[1])])
    # EPSG 4326 uses a coordinate system the same as a GLOBE (curved surface).
    # EPSG 3857 uses a coordinate system the same as a MAP (flat surface).
    project = partial(
        pyproj.transform,
        pyproj.Proj('EPSG:4326'),
        pyproj.Proj('EPSG:3857'))

    fixedLine = transform(project, line)
    return fixedLine.length


# print(distance([15.799406, 40.636069], [15.810173, 40.640246]))
# print(distance([15.799406, 40.636069], [15.810173, 40.640246]))