Convert bounding box values from decimal degrees to Web Mercator (Auxiliary Sphere) EPSG:3857

dd2webmercator.py

# This script converts bounding box values from decimal degrees to
# WGS 1984 Web Mercator (Auxiliary Sphere) EPSG:3857. Note that Web Mercator
# was originally assigned wkid 102100, but was later changed to 3857.
# See http://resources.arcgis.com/en/help/rest/apiref/geometry.html

# Arguments and Output
# With flags, arguments can be given in any order, but the output will always 
# be ordered like a GeoJSON bbox: "xmin, ymin, xmax, ymax". In both the GeoJSON and
# Esri JSON specs, a bounding box (or envelope) is defined with the lowest
# values for all axes followed by the highest values, e.g. "xmin, ymin, xmax, ymax".
# See http://geojson.org/geojson-spec.html#bounding-boxes
# See http://resources.arcgis.com/en/help/rest/apiref/geometry.html

# Example Usage
# dd2webmercator.py -w -94.292822041 -s 28.192227864 -e -88.645849385 -n 30.637075067
# dd2webmercator.py --west -94.292822041 --south 28.192227864 --east -88.645849385 --north 30.637075067

import math
import argparse

_rMajor = 6378137 # Equatorial Radius, WGS84
_shift = math.pi * _rMajor

def DecimalDegreesToWebMercator(lon, lat):
	x = lon * _shift / 180.0
	y = math.log(math.tan((90.0 + lat) * math.pi / 360.0)) / (math.pi / 180.0)
	y = y * _shift / 180.0
	return x,y
	
parser = argparse.ArgumentParser(description="Convert bounding box values from decimal degrees to 'Web Mercator (Auxiliary Sphere)' EPSG:3857", add_help=True)
group = parser.add_argument_group("Decimal Degrees Values")
group.add_argument("-w", "--west", "--west_longitude", "--xmin", type=float, required=True, help="West most longitude, aka 'minimum longitude', aka 'xmin', aka 'left'")
group.add_argument("-s", "--south", "--south_latitude", "--ymin", type=float, required=True, help="South most latitude, aka 'minimum latitude', aka 'ymin', aka 'bottom'")
group.add_argument("-e", "--east", "--east_longitude", "--xmax", type=float, required=True, help="East most longitude, aka 'maximum longitude', aka 'xmax', aka 'right'")
group.add_argument("-n", "--north", "--north_latitude", "--ymax", type=float, required=True, help="North most latitude, aka 'maximum latitude', aka 'ymax', aka 'top'")

args = parser.parse_args()

x1, y1 = DecimalDegreesToWebMercator(args.west, args.south)
x2, y2 = DecimalDegreesToWebMercator(args.east, args.north)

# Output values are: west, south, east, north
out = "%.10f,%.10f,%.10f,%.10f" % (x1, y1, x2, y2)

print out